110 class(
schwarz_t),
target,
intent(inout) :: this
111 type(
space_t),
target,
intent(inout) :: Xh
112 type(
dofmap_t),
target,
intent(in) :: dof
113 type(
gs_t),
target,
intent(inout) :: gs_h
114 type(
mesh_t),
target,
intent(inout) :: msh
115 type(
bc_list_t),
target,
intent(inout):: bclst
120 call this%Xh_schwarz%init(
gll, xh%lx+2, xh%lx+2, xh%lx+2)
121 call this%dm_schwarz%init(msh, this%Xh_schwarz)
122 call this%gs_schwarz%init(this%dm_schwarz)
124 allocate(this%work1(this%dm_schwarz%size()))
125 allocate(this%work2(this%dm_schwarz%size()))
126 allocate(this%wt(xh%lx, xh%lx, 4, msh%gdim, msh%nelv))
128 call this%fdm%init(xh, dof, gs_h)
143 if (nthrds .gt. 1)
then
145 call this%gs_h%init(this%dof)
146 this%local_gs = .true.
149 this%local_gs = .false.
153 call device_map(this%work1, this%work1_d,this%dm_schwarz%size())
154 call device_map(this%work2, this%work2_d,this%dm_schwarz%size())
160 int(this%dof%size() * c_sizeof(this%work1(1)),
i8))
161 call rone(this%work1, this%dof%size())
197 integer :: enx,eny,enz, n, ie, k, ns
198 real(kind=
rp),
parameter :: zero = 0.0
199 real(kind=
rp),
parameter :: one = 1.0
200 associate(work1 => this%work1, work2 => this%work2, msh => this%msh, &
201 xh => this%Xh, xh_schwarz => this%Xh_schwarz)
208 if(.not. msh%gdim .eq. 3) enz=1
209 ns = enx*eny*enz*msh%nelv
212 call rzero(work1, ns)
216 call schwarz_extrude(work1, 0, zero, work2, 0, one , enx, eny, enz, msh%nelv)
220 call this%gs_schwarz%op(work2, ns, gs_op_add)
224 call this%gs_schwarz%op(work2, ns, gs_op_add)
226 call schwarz_extrude(work2, 0, one, work1, 0, -one, enx, eny, enz, msh%nelv)
227 call schwarz_extrude(work2, 2, one, work2, 0, one, enx, eny, enz, msh%nelv)
238 call this%gs_h%op(work1, n, gs_op_add)
242 call this%gs_h%op(work1, n, gs_op_add)
247 if (msh%gdim .eq. 2)
then
250 if (this%msh%gdim.eq. 3)
then
260 integer,
intent(in) :: n, nelv
261 real(kind=rp),
intent(inout) :: wt(n,4,2,nelv)
262 real(kind=rp),
intent(inout) :: work(n,n)
265 wt(j,1,1,ie) = 1.0_rp / work(1,j)
266 wt(j,2,1,ie) = 1.0_rp / work(2,j)
267 wt(j,3,1,ie) = 1.0_rp / work(n-1,j)
268 wt(j,4,1,ie) = 1.0_rp / work(n,j)
271 wt(i,1,2,ie) = 1.0_rp / work(i,1)
272 wt(i,2,2,ie) = 1.0_rp / work(i,2)
273 wt(i,3,2,ie) = 1.0_rp / work(i,n-1)
274 wt(i,4,2,ie) = 1.0_rp / work(i,n)
282 integer,
intent(in) ::n, nelv, ie
283 real(kind=rp),
intent(inout) :: wt(n,n,4,3,nelv)
284 real(kind=rp),
intent(inout) :: work(n,n,n)
289 wt(j,k,1,1,ie) = 1.0_rp / work(1,j,k)
290 wt(j,k,2,1,ie) = 1.0_rp / work(2,j,k)
291 wt(j,k,3,1,ie) = 1.0_rp / work(n-1,j,k)
292 wt(j,k,4,1,ie) = 1.0_rp / work(n,j,k)
298 wt(i,k,1,2,ie) = 1.0_rp / work(i,1,k)
299 wt(i,k,2,2,ie) = 1.0_rp / work(i,2,k)
300 wt(i,k,3,2,ie) = 1.0_rp / work(i,n-1,k)
301 wt(i,k,4,2,ie) = 1.0_rp / work(i,n,k)
307 wt(i,j,1,3,ie) = 1.0_rp / work(i,j,1)
308 wt(i,j,2,3,ie) = 1.0_rp / work(i,j,2)
309 wt(i,j,3,3,ie) = 1.0_rp / work(i,j,n-1)
310 wt(i,j,4,3,ie) = 1.0_rp / work(i,j,n)
355 integer,
intent(in) :: l1, l2, nx, ny, nz, nelv
356 real(kind=rp),
intent(inout) :: arr1(nx,ny,nz,nelv), arr2(nx,ny,nz,nelv)
357 real(kind=rp),
intent(in) :: f1, f2
358 integer :: i, j, k, ie, i0, i1
365 arr1(l1+1 ,j,1,ie) = f1*arr1(l1+1 ,j,1,ie) &
366 +f2*arr2(l2+1 ,j,1,ie)
367 arr1(nx-l1,j,1,ie) = f1*arr1(nx-l1,j,1,ie) &
368 +f2*arr2(nx-l2,j,1,ie)
371 arr1(i,l1+1 ,1,ie) = f1*arr1(i,l1+1 ,1,ie) &
372 +f2*arr2(i,l2+1 ,1,ie)
373 arr1(i,ny-l1,1,ie) = f1*arr1(i,ny-l1,1,ie) &
374 +f2*arr2(i,nx-l2,1,ie)
381 arr1(l1+1 ,j,k,ie) = f1*arr1(l1+1 ,j,k,ie) &
382 +f2*arr2(l2+1 ,j,k,ie)
383 arr1(nx-l1,j,k,ie) = f1*arr1(nx-l1,j,k,ie) &
384 +f2*arr2(nx-l2,j,k,ie)
389 arr1(i,l1+1 ,k,ie) = f1*arr1(i,l1+1 ,k,ie) &
390 +f2*arr2(i,l2+1 ,k,ie)
391 arr1(i,nx-l1,k,ie) = f1*arr1(i,nx-l1,k,ie) &
392 +f2*arr2(i,nx-l2,k,ie)
397 arr1(i,j,l1+1 ,ie) = f1*arr1(i,j,l1+1 ,ie) &
398 +f2*arr2(i,j,l2+1 ,ie)
399 arr1(i,j,nx-l1,ie) = f1*arr1(i,j,nx-l1,ie) &
400 +f2*arr2(i,j,nx-l2,ie)
409 real(kind=rp),
dimension(this%dof%size()),
intent(inout) :: e, r
410 integer :: n, enx, eny, enz, ns
411 real(kind=rp),
parameter :: zero = 0.0_rp
412 real(kind=rp),
parameter :: one = 1.0_rp
413 type(c_ptr) :: e_d, r_d
414 associate(work1 => this%work1, work1_d => this%work1_d,&
415 work2 => this%work2, work2_d => this%work2_d)
418 enx=this%Xh_schwarz%lx
419 eny=this%Xh_schwarz%ly
420 enz=this%Xh_schwarz%lz
421 if(.not. this%msh%gdim .eq. 3) enz=1
422 ns = enx*eny*enz*this%msh%nelv
423 if (neko_bcknd_device .eq. 1)
then
424 r_d = device_get_ptr(r)
425 e_d = device_get_ptr(e)
426 call device_event_record(this%event, glb_cmd_queue)
427 call device_stream_wait_event(aux_cmd_queue, this%event, 0)
428 call device_schwarz_toext3d(work1_d, r_d, this%Xh%lx, &
429 this%msh%nelv, aux_cmd_queue)
430 call device_schwarz_extrude(work1_d, 0, zero, work1_d, 2, one, &
431 enx,eny,enz, this%msh%nelv,aux_cmd_queue)
433 this%gs_schwarz%bcknd%gs_stream = aux_cmd_queue
434 call this%gs_schwarz%op(work1, ns, gs_op_add,this%event)
435 call device_event_sync(this%event)
436 call device_schwarz_extrude(work1_d, 0, one, work1_d, 2, -one, &
437 enx, eny, enz, this%msh%nelv, aux_cmd_queue)
439 call this%fdm%compute(work2, work1,aux_cmd_queue)
441 call device_schwarz_extrude(work1_d, 0, zero, work2_d, 0, one, &
442 enx, eny, enz, this%msh%nelv, aux_cmd_queue)
443 call this%gs_schwarz%op(work2, ns, gs_op_add,this%event)
444 call device_event_sync(this%event)
446 call device_schwarz_extrude(work2_d, 0, one, work1_d, 0, -one, &
447 enx, eny, enz, this%msh%nelv, aux_cmd_queue)
448 call device_schwarz_extrude(work2_d, 2, one, work2_d, 0, one, &
449 enx, eny, enz, this%msh%nelv, aux_cmd_queue)
450 call device_schwarz_toreg3d(e_d, work2_d, this%Xh%lx, &
451 this%msh%nelv, aux_cmd_queue)
453 this%gs_h%bcknd%gs_stream = aux_cmd_queue
454 call this%gs_h%op(e, n, gs_op_add, this%event)
456 call this%bclst%apply_scalar(e, n, strm = aux_cmd_queue)
457 call device_col2(e_d,this%wt_d, n, aux_cmd_queue)
460 if (.not. this%local_gs)
then
461 call device_event_sync(this%event)
462 this%gs_h%bcknd%gs_stream = glb_cmd_queue
465 call this%bclst%apply_scalar(r, n)
470 enx, eny, enz, this%msh%nelv)
471 call this%gs_schwarz%op(work1, ns, gs_op_add)
473 enx, eny, enz, this%msh%nelv)
475 call this%fdm%compute(work2, work1)
479 enx, eny, enz, this%msh%nelv)
480 call this%gs_schwarz%op(work2, ns, gs_op_add)
482 enx, eny, enz, this%msh%nelv)
484 enx, eny, enz, this%msh%nelv)
489 call this%gs_h%op(e, n, gs_op_add)
490 call this%bclst%apply_scalar(e, n)
492 call schwarz_wt3d(e, this%wt, this%Xh%lx, this%msh%nelv)
499 integer,
intent(in) :: n, nelv
500 real(kind=rp),
intent(inout) :: e(n,n,n,nelv)
501 real(kind=rp),
intent(inout) :: wt(n,n,4,3,nelv)
502 integer :: ie, i, j, k
507 e(1,j,k,ie) = e(1,j,k,ie) * wt(j,k,1,1,ie)
508 e(2,j,k,ie) = e(2,j,k,ie) * wt(j,k,2,1,ie)
509 e(n-1,j,k,ie) = e(n-1,j,k,ie) * wt(j,k,3,1,ie)
510 e(n,j,k,ie) = e(n,j,k,ie) * wt(j,k,4,1,ie)
515 e(i,1,k,ie) = e(i,1,k,ie) * wt(i,k,1,2,ie)
516 e(i,2,k,ie) = e(i,2,k,ie) * wt(i,k,2,2,ie)
517 e(i,n-1,k,ie) = e(i,n-1,k,ie) * wt(i,k,3,2,ie)
518 e(i,n,k,ie) = e(i,n,k,ie) * wt(i,k,4,2,ie)
523 e(i,j,1,ie) = e(i,j,1,ie) * wt(i,j,1,3,ie)
524 e(i,j,2,ie) = e(i,j,2,ie) * wt(i,j,2,3,ie)
525 e(i,j,n-1,ie) = e(i,j,n-1,ie) * wt(i,j,3,3,ie)
526 e(i,j,n,ie) = e(i,j,n,ie) * wt(i,j,4,3,ie)