Neko 0.9.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
schwarz.f90
Go to the documentation of this file.
1! Copyright (c) 2008-2020, UCHICAGO ARGONNE, LLC.
2!
3! The UChicago Argonne, LLC as Operator of Argonne National
4! Laboratory holds copyright in the Software. The copyright holder
5! reserves all rights except those expressly granted to licensees,
6! and U.S. Government license rights.
7!
8! Redistribution and use in source and binary forms, with or without
9! modification, are permitted provided that the following conditions
10! are met:
11!
12! 1. Redistributions of source code must retain the above copyright
13! notice, this list of conditions and the disclaimer below.
14!
15! 2. Redistributions in binary form must reproduce the above copyright
16! notice, this list of conditions and the disclaimer (as noted below)
17! in the documentation and/or other materials provided with the
18! distribution.
19!
20! 3. Neither the name of ANL nor the names of its contributors
21! may be used to endorse or promote products derived from this software
22! without specific prior written permission.
23!
24! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
25! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
26! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
27! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL
28! UCHICAGO ARGONNE, LLC, THE U.S. DEPARTMENT OF
29! ENERGY OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
30! SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
31! TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
32! DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
33! THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
34! (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
35! OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36!
37! Additional BSD Notice
38! ---------------------
39! 1. This notice is required to be provided under our contract with
40! the U.S. Department of Energy (DOE). This work was produced at
41! Argonne National Laboratory under Contract
42! No. DE-AC02-06CH11357 with the DOE.
43!
44! 2. Neither the United States Government nor UCHICAGO ARGONNE,
45! LLC nor any of their employees, makes any warranty,
46! express or implied, or assumes any liability or responsibility for the
47! accuracy, completeness, or usefulness of any information, apparatus,
48! product, or process disclosed, or represents that its use would not
49! infringe privately-owned rights.
50!
51! 3. Also, reference herein to any specific commercial products, process,
52! or services by trade name, trademark, manufacturer or otherwise does
53! not necessarily constitute or imply its endorsement, recommendation,
54! or favoring by the United States Government or UCHICAGO ARGONNE LLC.
55! The views and opinions of authors expressed
56! herein do not necessarily state or reflect those of the United States
57! Government or UCHICAGO ARGONNE, LLC, and shall
58! not be used for advertising or product endorsement purposes.
59!
61module schwarz
62 use num_types
63 use math
64 use mesh
65 use space
66 use dofmap, only : dofmap_t
67 use bc
70 use device_math
71 use fdm, only : fdm_t
72 use device
73 use neko_config
74 use, intrinsic :: iso_c_binding
75 implicit none
76 private
77
78 type, public :: schwarz_t
79 real(kind=rp), allocatable :: work1(:)
80 real(kind=rp), allocatable :: work2(:)
81 real(kind=rp), allocatable :: wt(:,:,:,:,:)
82 type(c_ptr) :: work1_d = c_null_ptr
83 type(c_ptr) :: work2_d = c_null_ptr
84 type(c_ptr) :: wt_d = c_null_ptr
85 type(space_t) :: xh_schwarz
86 type(gs_t) :: gs_schwarz
87 type(dofmap_t) :: dm_schwarz
88 type(fdm_t) :: fdm
89 type(space_t), pointer :: xh
90 type(bc_list_t), pointer :: bclst
91 type(dofmap_t), pointer :: dm
92 type(gs_t), pointer :: gs_h
93 type(mesh_t), pointer :: msh
94 type(c_ptr) :: event
95 contains
96 procedure, pass(this) :: init => schwarz_init
97 procedure, pass(this) :: free => schwarz_free
98 procedure, pass(this) :: compute => schwarz_compute
99 end type schwarz_t
100
101contains
102
103 subroutine schwarz_init(this, Xh, dm, gs_h, bclst, msh)
104 class(schwarz_t), target, intent(inout) :: this
105 type(space_t), target, intent(inout) :: Xh
106 type(dofmap_t), target, intent(inout) :: dm
107 type(gs_t), target, intent(inout) :: gs_h
108 type(mesh_t), target, intent(inout) :: msh
109 type(bc_list_t), target, intent(inout):: bclst
110
111 call this%free()
112
113 call this%Xh_schwarz%init(gll, xh%lx+2, xh%lx+2, xh%lx+2)
114 call this%dm_schwarz%init(msh, this%Xh_schwarz)
115 call this%gs_schwarz%init(this%dm_schwarz)
116
117 allocate(this%work1(this%dm_schwarz%size()))
118 allocate(this%work2(this%dm_schwarz%size()))
119 allocate(this%wt(xh%lx, xh%lx, 4, msh%gdim, msh%nelv))
120
121 call this%fdm%init(xh, dm, gs_h)
122
123
124 this%msh => msh
125 this%Xh => xh
126 this%bclst => bclst
127 this%dm => dm
128 this%gs_h => gs_h
129 if (neko_bcknd_device .eq. 1) then
130 call device_map(this%work1, this%work1_d,this%dm_schwarz%size())
131 call device_map(this%work2, this%work2_d,this%dm_schwarz%size())
132 end if
133
134
135 call schwarz_setup_wt(this)
136 if (neko_bcknd_device .eq. 1) then
137 call device_alloc(this%wt_d,int(this%dm%size()*rp, i8))
138 call rone(this%work1, this%dm%size())
139 call schwarz_wt3d(this%work1, this%wt, xh%lx, msh%nelv)
140 call device_memcpy(this%work1, this%wt_d, this%dm%size(), &
141 host_to_device, sync=.false.)
142 call device_event_create(this%event, 2)
143 end if
144 end subroutine schwarz_init
145
146 subroutine schwarz_free(this)
147 class(schwarz_t), intent(inout) :: this
148
149 if(allocated(this%work1)) deallocate(this%work1)
150 if(allocated(this%work2)) deallocate(this%work2)
151 if(allocated(this%wt)) deallocate(this%wt)
152
153 call this%Xh_schwarz%free()
154 call this%gs_schwarz%free()
155 !why cant I do this?
156 !call dofmap_free(this%dm_schwarz)
157 call this%fdm%free()
158
159 nullify(this%Xh)
160 nullify(this%bclst)
161 nullify(this%dm)
162 nullify(this%gs_h)
163 nullify(this%msh)
164 end subroutine schwarz_free
166 subroutine schwarz_setup_wt(this)
167 class(schwarz_t), intent(inout) :: this
168 integer :: enx,eny,enz, n, ie, k, ns
169 real(kind=rp), parameter :: zero = 0.0
170 real(kind=rp), parameter :: one = 1.0
171 associate(work1 => this%work1, work2 => this%work2, msh => this%msh, &
172 xh => this%Xh, xh_schwarz => this%Xh_schwarz)
173
174 n = this%dm%size()
175
176 enx = xh_schwarz%lx
177 eny = xh_schwarz%ly
178 enz = xh_schwarz%lz
179 if(.not. msh%gdim .eq. 3) enz=1
180 ns = enx*eny*enz*msh%nelv
181
182 call rone(work2, ns)
183 call rzero(work1, ns)
184
185 ! Sum overlap region (border excluded)
186 ! Cred to PFF for this, very clever
187 call schwarz_extrude(work1, 0, zero, work2, 0, one , enx, eny, enz, msh%nelv)
188 if (neko_bcknd_device .eq. 1) then
189 call device_memcpy(work2, this%work2_d, ns, &
190 host_to_device, sync=.false.)
191 call this%gs_schwarz%op(work2, ns, gs_op_add)
192 call device_memcpy(work2, this%work2_d, ns, &
193 device_to_host, sync=.false.)
194 else
195 call this%gs_schwarz%op(work2, ns, gs_op_add)
196 end if
197 call schwarz_extrude(work2, 0, one, work1, 0, -one, enx, eny, enz, msh%nelv)
198 call schwarz_extrude(work2, 2, one, work2, 0, one, enx, eny, enz, msh%nelv)
199
200 ! if(.not.if3d) then ! Go back to regular size array
201 ! call hsmg_schwarz_toreg2d(mg_work,mg_work(i),mg_nh(l))
202 ! else
203 call schwarz_toreg3d(work1, work2, xh%lx, msh%nelv)
204 ! endif
205
206 if (neko_bcknd_device .eq. 1) then
207 call device_memcpy(work1, this%work1_d, n, &
208 host_to_device, sync=.false.)
209 call this%gs_h%op(work1, n, gs_op_add)
210 call device_memcpy(work1, this%work1_d, n, &
211 device_to_host, sync=.true.)
212 else
213 call this%gs_h%op(work1, n, gs_op_add)
214 end if
215
216 k = 1
217 do ie = 1,msh%nelv
218 if (msh%gdim .eq. 2) then
219 call schwarz_setup_schwarz_wt2d_2(this%wt,ie,xh%lx, work1(k), msh%nelv)
220 end if
221 if (this%msh%gdim.eq. 3) then
222 call schwarz_setup_schwarz_wt3d_2(this%wt,ie,xh%lx, work1(k), msh%nelv)
223 k = k + xh%lxyz
224 end if
225 end do
226 end associate
227 end subroutine schwarz_setup_wt
228
230 subroutine schwarz_setup_schwarz_wt2d_2(wt,ie,n,work, nelv)
231 integer, intent(in) :: n, nelv
232 real(kind=rp), intent(inout) :: wt(n,4,2,nelv)
233 real(kind=rp), intent(inout) :: work(n,n)
234 integer :: ie,i,j
235 do j = 1, n
236 wt(j,1,1,ie) = 1.0_rp / work(1,j)
237 wt(j,2,1,ie) = 1.0_rp / work(2,j)
238 wt(j,3,1,ie) = 1.0_rp / work(n-1,j)
239 wt(j,4,1,ie) = 1.0_rp / work(n,j)
240 end do
241 do i = 1, n
242 wt(i,1,2,ie) = 1.0_rp / work(i,1)
243 wt(i,2,2,ie) = 1.0_rp / work(i,2)
244 wt(i,3,2,ie) = 1.0_rp / work(i,n-1)
245 wt(i,4,2,ie) = 1.0_rp / work(i,n)
246 end do
247
248 return
249 end subroutine schwarz_setup_schwarz_wt2d_2
250
252 subroutine schwarz_setup_schwarz_wt3d_2(wt, ie, n, work, nelv)
253 integer, intent(in) ::n, nelv, ie
254 real(kind=rp), intent(inout) :: wt(n,n,4,3,nelv)
255 real(kind=rp), intent(inout) :: work(n,n,n)
256 integer :: i,j,k
257
258 do k = 1, n
259 do j = 1, n
260 wt(j,k,1,1,ie) = 1.0_rp / work(1,j,k)
261 wt(j,k,2,1,ie) = 1.0_rp / work(2,j,k)
262 wt(j,k,3,1,ie) = 1.0_rp / work(n-1,j,k)
263 wt(j,k,4,1,ie) = 1.0_rp / work(n,j,k)
264 end do
265 end do
266
267 do k = 1, n
268 do i = 1, n
269 wt(i,k,1,2,ie) = 1.0_rp / work(i,1,k)
270 wt(i,k,2,2,ie) = 1.0_rp / work(i,2,k)
271 wt(i,k,3,2,ie) = 1.0_rp / work(i,n-1,k)
272 wt(i,k,4,2,ie) = 1.0_rp / work(i,n,k)
273 end do
274 end do
275
276 do j = 1, n
277 do i = 1, n
278 wt(i,j,1,3,ie) = 1.0_rp / work(i,j,1)
279 wt(i,j,2,3,ie) = 1.0_rp / work(i,j,2)
280 wt(i,j,3,3,ie) = 1.0_rp / work(i,j,n-1)
281 wt(i,j,4,3,ie) = 1.0_rp / work(i,j,n)
282 end do
283 end do
284
285 end subroutine schwarz_setup_schwarz_wt3d_2
286
288 subroutine schwarz_toreg3d(b, a, n, nelv)
289 integer, intent(in) :: n, nelv
290 real(kind=rp), intent(inout) :: a(0:n+1, 0:n+1, 0:n+1, nelv)
291 real(kind=rp), intent(inout) :: b(n,n,n,nelv)
292 integer :: i, j, k, ie
293 do ie = 1, nelv
294 do k = 1, n
295 do j = 1, n
296 do i = 1, n
297 b(i,j,k,ie) = a(i,j,k,ie)
298 end do
299 end do
300 end do
301 end do
302 end subroutine schwarz_toreg3d
303
305 subroutine schwarz_toext3d(a, b, n, nelv)
306 integer, intent(in) :: n, nelv
307 real (kind=rp), intent(inout) :: a(0:n+1,0:n+1,0:n+1,nelv), b(n,n,n,nelv)
308 integer :: i,j,k,ie
309
310 call rzero(a, (n+2)*(n+2)*(n+2)*nelv)
311 do ie = 1, nelv
312 do k = 1, n
313 do j = 1, n
314 do i = 1, n
315 a(i,j,k,ie) = b(i,j,k,ie)
316 end do
317 end do
318 end do
319 end do
320 end subroutine schwarz_toext3d
321
325 subroutine schwarz_extrude(arr1, l1, f1, arr2, l2, f2, nx, ny, nz, nelv)
326 integer, intent(in) :: l1, l2, nx, ny, nz, nelv
327 real(kind=rp), intent(inout) :: arr1(nx,ny,nz,nelv), arr2(nx,ny,nz,nelv)
328 real(kind=rp), intent(in) :: f1, f2
329 integer :: i, j, k, ie, i0, i1
330 i0=2
331 i1=nx-1
332
333 if(nz .eq. 1) then
334 do ie = 1, nelv
335 do j = i0, i1
336 arr1(l1+1 ,j,1,ie) = f1*arr1(l1+1 ,j,1,ie) &
337 +f2*arr2(l2+1 ,j,1,ie)
338 arr1(nx-l1,j,1,ie) = f1*arr1(nx-l1,j,1,ie) &
339 +f2*arr2(nx-l2,j,1,ie)
340 end do
341 do i = i0, i1
342 arr1(i,l1+1 ,1,ie) = f1*arr1(i,l1+1 ,1,ie) &
343 +f2*arr2(i,l2+1 ,1,ie)
344 arr1(i,ny-l1,1,ie) = f1*arr1(i,ny-l1,1,ie) &
345 +f2*arr2(i,nx-l2,1,ie)
346 end do
347 end do
348 else
349 do ie = 1, nelv
350 do k = i0, i1
351 do j = i0, i1
352 arr1(l1+1 ,j,k,ie) = f1*arr1(l1+1 ,j,k,ie) &
353 +f2*arr2(l2+1 ,j,k,ie)
354 arr1(nx-l1,j,k,ie) = f1*arr1(nx-l1,j,k,ie) &
355 +f2*arr2(nx-l2,j,k,ie)
356 end do
357 end do
358 do k = i0, i1
359 do i = i0, i1
360 arr1(i,l1+1 ,k,ie) = f1*arr1(i,l1+1 ,k,ie) &
361 +f2*arr2(i,l2+1 ,k,ie)
362 arr1(i,nx-l1,k,ie) = f1*arr1(i,nx-l1,k,ie) &
363 +f2*arr2(i,nx-l2,k,ie)
364 end do
365 end do
366 do j = i0, i1
367 do i = i0, i1
368 arr1(i,j,l1+1 ,ie) = f1*arr1(i,j,l1+1 ,ie) &
369 +f2*arr2(i,j,l2+1 ,ie)
370 arr1(i,j,nx-l1,ie) = f1*arr1(i,j,nx-l1,ie) &
371 +f2*arr2(i,j,nx-l2,ie)
372 end do
373 end do
374 end do
375 endif
376 end subroutine schwarz_extrude
377
378 subroutine schwarz_compute(this, e, r)
379 class(schwarz_t), intent(inout) :: this
380 real(kind=rp), dimension(this%dm%size()), intent(inout) :: e, r
381 integer :: n, enx, eny, enz, ns
382 real(kind=rp), parameter :: zero = 0.0_rp
383 real(kind=rp), parameter :: one = 1.0_rp
384 type(c_ptr) :: e_d, r_d
385 associate(work1 => this%work1, work1_d => this%work1_d,&
386 work2 => this%work2, work2_d => this%work2_d)
387
388 n = this%dm%size()
389 enx=this%Xh_schwarz%lx
390 eny=this%Xh_schwarz%ly
391 enz=this%Xh_schwarz%lz
392 if(.not. this%msh%gdim .eq. 3) enz=1
393 ns = enx*eny*enz*this%msh%nelv
394 if (neko_bcknd_device .eq. 1) then
395 r_d = device_get_ptr(r)
396 e_d = device_get_ptr(e)
397 call device_event_record(this%event, glb_cmd_queue)
398 call device_stream_wait_event(aux_cmd_queue, this%event, 0)
399 call device_schwarz_toext3d(work1_d, r_d, this%Xh%lx, &
400 this%msh%nelv, aux_cmd_queue)
401 call device_schwarz_extrude(work1_d, 0, zero, work1_d, 2, one, &
402 enx,eny,enz, this%msh%nelv,aux_cmd_queue)
403
404 this%gs_schwarz%bcknd%gs_stream = aux_cmd_queue
405 call this%gs_schwarz%op(work1, ns, gs_op_add,this%event)
406 call device_event_sync(this%event)
407 call device_schwarz_extrude(work1_d, 0, one, work1_d, 2, -one, &
408 enx, eny, enz, this%msh%nelv, aux_cmd_queue)
409
410 call this%fdm%compute(work2, work1,aux_cmd_queue) ! do local solves
411
412 call device_schwarz_extrude(work1_d, 0, zero, work2_d, 0, one, &
413 enx, eny, enz, this%msh%nelv, aux_cmd_queue)
414 call this%gs_schwarz%op(work2, ns, gs_op_add,this%event)
415 call device_event_sync(this%event)
416
417 call device_schwarz_extrude(work2_d, 0, one, work1_d, 0, -one, &
418 enx, eny, enz, this%msh%nelv, aux_cmd_queue)
419 call device_schwarz_extrude(work2_d, 2, one, work2_d, 0, one, &
420 enx, eny, enz, this%msh%nelv, aux_cmd_queue)
421 call device_schwarz_toreg3d(e_d, work2_d, this%Xh%lx, &
422 this%msh%nelv, aux_cmd_queue)
423
424 call device_event_record(this%event,aux_cmd_queue)
425 call device_event_sync(this%event)
426
427 call this%gs_h%op(e, n, gs_op_add, this%event)
428 call bc_list_apply_scalar(this%bclst, e, n)
429 call device_col2(e_d,this%wt_d, n)
430 call device_stream_wait_event(aux_cmd_queue, this%event, 0)
431 else
432 call bc_list_apply_scalar(this%bclst, r, n)
433 call schwarz_toext3d(work1, r, this%Xh%lx, this%msh%nelv)
434
435 ! exchange interior nodes
436 call schwarz_extrude(work1, 0, zero, work1, 2, one, &
437 enx, eny, enz, this%msh%nelv)
438 call this%gs_schwarz%op(work1, ns, gs_op_add)
439 call schwarz_extrude(work1, 0, one, work1, 2, -one, &
440 enx, eny, enz, this%msh%nelv)
441
442 call this%fdm%compute(work2, work1) ! do local solves
443
444 ! Sum overlap region (border excluded)
445 call schwarz_extrude(work1, 0, zero, work2, 0, one, &
446 enx, eny, enz, this%msh%nelv)
447 call this%gs_schwarz%op(work2, ns, gs_op_add)
448 call schwarz_extrude(work2, 0, one, work1, 0, -one, &
449 enx, eny, enz, this%msh%nelv)
450 call schwarz_extrude(work2, 2, one, work2, 0, one, &
451 enx, eny, enz, this%msh%nelv)
452
453 call schwarz_toreg3d(e, work2, this%Xh%lx, this%msh%nelv)
454
455 ! sum border nodes
456 call this%gs_h%op(e, n, gs_op_add)
457 call bc_list_apply_scalar(this%bclst, e, n)
458
459 call schwarz_wt3d(e, this%wt, this%Xh%lx, this%msh%nelv)
460 end if
461 end associate
462 end subroutine schwarz_compute
463
464 !Apply schwarz weights along the boundary of each element.
465 subroutine schwarz_wt3d(e,wt,n, nelv)
466 integer, intent(in) :: n, nelv
467 real(kind=rp), intent(inout) :: e(n,n,n,nelv)
468 real(kind=rp), intent(inout) :: wt(n,n,4,3,nelv)
469 integer :: ie, i, j, k
470
471 do ie = 1, nelv
472 do k = 1, n
473 do j = 1, n
474 e(1 ,j,k,ie) = e(1 ,j,k,ie) * wt(j,k,1,1,ie)
475 e(2 ,j,k,ie) = e(2 ,j,k,ie) * wt(j,k,2,1,ie)
476 e(n-1,j,k,ie) = e(n-1,j,k,ie) * wt(j,k,3,1,ie)
477 e(n ,j,k,ie) = e(n ,j,k,ie) * wt(j,k,4,1,ie)
478 end do
479 end do
480 do k = 1, n
481 do i = 3, n-2
482 e(i,1 ,k,ie) = e(i,1 ,k,ie) * wt(i,k,1,2,ie)
483 e(i,2 ,k,ie) = e(i,2 ,k,ie) * wt(i,k,2,2,ie)
484 e(i,n-1,k,ie) = e(i,n-1,k,ie) * wt(i,k,3,2,ie)
485 e(i,n ,k,ie) = e(i,n ,k,ie) * wt(i,k,4,2,ie)
486 end do
487 end do
488 do j = 3, n-2
489 do i = 3, n-2
490 e(i,j,1 ,ie) = e(i,j,1 ,ie) * wt(i,j,1,3,ie)
491 e(i,j,2 ,ie) = e(i,j,2 ,ie) * wt(i,j,2,3,ie)
492 e(i,j,n-1,ie) = e(i,j,n-1,ie) * wt(i,j,3,3,ie)
493 e(i,j,n ,ie) = e(i,j,n ,ie) * wt(i,j,4,3,ie)
494 end do
495 end do
496 end do
497 end subroutine schwarz_wt3d
498end module schwarz
Map a Fortran array to a device (allocate and associate)
Definition device.F90:57
Copy data between host and device (or device and device)
Definition device.F90:51
Defines a boundary condition.
Definition bc.f90:34
Device abstraction, common interface for various accelerators.
Definition device.F90:34
integer, parameter, public host_to_device
Definition device.F90:47
integer, parameter, public device_to_host
Definition device.F90:47
subroutine, public device_alloc(x_d, s)
Allocate memory on the device.
Definition device.F90:164
subroutine, public device_event_create(event, flags)
Create a device event queue.
Definition device.F90:1164
Defines a mapping of the degrees of freedom.
Definition dofmap.f90:35
Type for the Fast Diagonalization connected with the schwarz overlapping solves.
Definition fdm.f90:61
Gather-scatter.
Definition math.f90:60
subroutine, public rone(a, n)
Set all elements to one.
Definition math.f90:228
subroutine, public rzero(a, n)
Zero a real vector.
Definition math.f90:195
Defines a mesh.
Definition mesh.f90:34
Build configurations.
integer, parameter neko_bcknd_device
integer, parameter, public i8
Definition num_types.f90:7
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Overlapping schwarz solves.
Definition schwarz.f90:61
subroutine schwarz_toreg3d(b, a, n, nelv)
convert array a from extended size to regular
Definition schwarz.f90:289
subroutine schwarz_setup_wt(this)
setup weights
Definition schwarz.f90:167
subroutine schwarz_extrude(arr1, l1, f1, arr2, l2, f2, nx, ny, nz, nelv)
Sum values along rows l1, l2 with weights f1, f2 and store along row l1. Helps us avoid complicated c...
Definition schwarz.f90:326
subroutine schwarz_setup_schwarz_wt3d_2(wt, ie, n, work, nelv)
Setup schwarz weights, 3d, second step.
Definition schwarz.f90:253
subroutine schwarz_compute(this, e, r)
Definition schwarz.f90:379
subroutine schwarz_wt3d(e, wt, n, nelv)
Definition schwarz.f90:466
subroutine schwarz_toext3d(a, b, n, nelv)
convert array a from original size to size extended array with border
Definition schwarz.f90:306
subroutine schwarz_init(this, xh, dm, gs_h, bclst, msh)
Definition schwarz.f90:104
subroutine schwarz_setup_schwarz_wt2d_2(wt, ie, n, work, nelv)
Setup schwarz weights, 2d, second step.
Definition schwarz.f90:231
subroutine schwarz_free(this)
Definition schwarz.f90:147
Defines a function space.
Definition space.f90:34
integer, parameter, public gll
Definition space.f90:48
A list of boundary conditions.
Definition bc.f90:104
The function space for the SEM solution fields.
Definition space.f90:62