Fills the extBoundary attribute on all faces participating in a 2:1 nonconforming (mortar) interface; vector analogue of the scalar MortarExchange (see MappedScalar3D_t for the algorithm description).
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| class(MappedVector3D_t), | intent(inout) | :: | this | |||
| type(Mesh3D), | intent(inout) | :: | mesh |
subroutine MortarExchange_MappedVector3D_t(this,mesh)
!! Fills the extBoundary attribute on all faces participating in a 2:1
!! nonconforming (mortar) interface; vector analogue of the scalar MortarExchange
!! (see MappedScalar3D_t for the algorithm description).
implicit none
class(MappedVector3D_t),intent(inout) :: this
type(Mesh3D),intent(inout) :: mesh
! Local
integer :: m,q,ivar,idir,i,j,i2,j2
integer :: eB,eS,sS,flip
integer :: rankId,offset,N
integer,pointer :: elemtorank(:)
real(prec) :: extBuff(1:this%interp%N+1,1:this%interp%N+1)
! See https://github.com/FluidNumerics/SELF/issues/54 for the reason behind
! this pointer alias
elemtorank => mesh%decomp%elemToRank(:)
rankId = mesh%decomp%rankId
offset = mesh%decomp%offsetElem(rankId+1)
N = this%interp%N
if(.not. allocated(this%mortarBuff)) then
allocate(this%mortarBuff(1:N+1,1:N+1,1:8,1:mesh%nMortars,1:this%nvar,1:3))
this%mortarBuff = 0.0_prec
endif
if(mesh%decomp%mpiEnabled) then
call this%MPIMortarExchangeAsync(mesh)
endif
! Stage rank-local traces in the big face's coordinates
do concurrent(m=1:mesh%nMortars,ivar=1:this%nvar,idir=1:3)
block
integer :: i,j,i2,j2,q
integer :: eB,sB,eS,sS,flip
eB = mesh%mortarInfo(1,m)
if(elemtorank(eB) == rankId) then
sB = mesh%mortarInfo(2,m)
do q = 1,4
do j = 1,N+1
do i = 1,N+1
this%mortarBuff(i,j,q,m,ivar,idir) = &
this%boundary(i,j,sB,eB-offset,ivar,idir)
enddo
enddo
enddo
endif
do q = 1,4
eS = mesh%mortarInfo(2*q+1,m)
if(elemtorank(eS) == rankId) then
sS = mesh%mortarInfo(2*q+2,m)/10
flip = mesh%mortarInfo(2*q+2,m)-10*sS
do j = 1,N+1
do i = 1,N+1
call MortarFaceMap(i,j,N,flip,i2,j2)
this%mortarBuff(i,j,4+q,m,ivar,idir) = &
this%boundary(i2,j2,sS,eS-offset,ivar,idir)
enddo
enddo
endif
enddo
endblock
enddo
if(mesh%decomp%mpiEnabled) then
call mesh%decomp%FinalizeMPIExchangeAsync()
! Reorient small-face traces received over MPI into the big face's coordinates
do idir = 1,3
do ivar = 1,this%nvar
do m = 1,mesh%nMortars
eB = mesh%mortarInfo(1,m)
if(elemtorank(eB) == rankId) then
do q = 1,4
eS = mesh%mortarInfo(2*q+1,m)
sS = mesh%mortarInfo(2*q+2,m)/10
flip = mesh%mortarInfo(2*q+2,m)-10*sS
if(elemtorank(eS) /= rankId .and. flip /= 0) then
do j = 1,N+1
do i = 1,N+1
call MortarFaceMap(i,j,N,flip,i2,j2)
extBuff(i,j) = this%mortarBuff(i2,j2,4+q,m,ivar,idir)
enddo
enddo
do j = 1,N+1
do i = 1,N+1
this%mortarBuff(i,j,4+q,m,ivar,idir) = extBuff(i,j)
enddo
enddo
endif
enddo
endif
enddo
enddo
enddo
endif
! Compute external states :
! small faces get the restricted big-face trace (exact),
! the big face gets the L2 projection of the small-face traces
do concurrent(m=1:mesh%nMortars,ivar=1:this%nvar,idir=1:3)
block
integer :: i,j,ii,jj,i2,j2,kx,ky,q
integer :: eB,sB,eS,sS,flip
real(prec) :: fm
real(prec) :: tmp(1:N+1,1:N+1)
real(prec) :: acc(1:N+1,1:N+1)
do q = 1,4
eS = mesh%mortarInfo(2*q+1,m)
if(elemtorank(eS) == rankId) then
sS = mesh%mortarInfo(2*q+2,m)/10
flip = mesh%mortarInfo(2*q+2,m)-10*sS
kx = mortarQuadKx(q)
ky = mortarQuadKy(q)
do jj = 1,N+1
do i = 1,N+1
fm = 0.0_prec
do ii = 1,N+1
fm = fm+this%interp%mortarR(ii,i,kx)* &
this%mortarBuff(ii,jj,q,m,ivar,idir)
enddo
tmp(i,jj) = fm
enddo
enddo
do j = 1,N+1
do i = 1,N+1
fm = 0.0_prec
do jj = 1,N+1
fm = fm+this%interp%mortarR(jj,j,ky)*tmp(i,jj)
enddo
call MortarFaceMap(i,j,N,flip,i2,j2)
this%extBoundary(i2,j2,sS,eS-offset,ivar,idir) = fm
enddo
enddo
endif
enddo
eB = mesh%mortarInfo(1,m)
if(elemtorank(eB) == rankId) then
sB = mesh%mortarInfo(2,m)
acc = 0.0_prec
do q = 1,4
kx = mortarQuadKx(q)
ky = mortarQuadKy(q)
do jj = 1,N+1
do i = 1,N+1
fm = 0.0_prec
do ii = 1,N+1
fm = fm+this%interp%mortarP(ii,i,kx)* &
this%mortarBuff(ii,jj,4+q,m,ivar,idir)
enddo
tmp(i,jj) = fm
enddo
enddo
do j = 1,N+1
do i = 1,N+1
fm = 0.0_prec
do jj = 1,N+1
fm = fm+this%interp%mortarP(jj,j,ky)*tmp(i,jj)
enddo
acc(i,j) = acc(i,j)+fm
enddo
enddo
enddo
do j = 1,N+1
do i = 1,N+1
this%extBoundary(i,j,sB,eB-offset,ivar,idir) = acc(i,j)
enddo
enddo
endif
endblock
enddo
endsubroutine MortarExchange_MappedVector3D_t