1 subroutine qpbc( rmass, qres, iBC, iper, ilwork ) 2c--------------------------------------------------------------------- 3c 4c This routine satisfies the periodic boundary conditions 5c on the diffusive flux residual and mass matrix 6c 7c input: 8c rmass (nshg) : mass matrix 9c qres (nshg,(nflow-1)*nsd) : diffusive flux vector 10c 11c output: modified qres and rmass 12c--------------------------------------------------------------------- 13 include "common.h" 14 15 dimension rmass(nshg), qres(nshg,idflx), 16 & iBC(nshg), iper(nshg),uv(nshg,2), 17 & tmpvec(nshg,4), tmp(nshg) 18c 19 if(iabc==1) then !are there any axisym bc's 20 do i=1,idflx/nsd 21 do j=1,2 22 istrt=j+(i-1)*(nflow-1) 23 uv(:,j)=qres(:,istrt) 24 enddo 25 call rotabc(uv, iBC, 'in ') 26 do j=1,2 27 istrt=j+(i-1)*(nflow-1) 28 qres(:,istrt)=uv(:,j) 29 enddo 30 enddo 31 endif 32c 33c 34c.... compute qi for node A, i.e., qres <-- qres/rmass 35c 36 if (numpe > 1) then 37 call commu (qres , ilwork, idflx , 'in ') 38 call commu (rmass , ilwork, 1 , 'in ') 39 endif 40c 41c take care of periodic boundary conditions 42c but not on surface tension terms in qres(:,10-12) 43c that are used to compute normal vector 44c 45 idflow = (nflow-1)*nsd 46 do j= 1,nshg 47 if ((btest(iBC(j),10))) then 48 i = iper(j) 49 rmass(i) = rmass(i) + rmass(j) 50c qres(i,:) = qres(i,:) + qres(j,:) 51 qres(i,1:idflow) = qres(i,1:idflow) + qres(j,1:idflow) 52 endif 53 enddo 54 55 do j= 1,nshg 56 if ((btest(iBC(j),10))) then 57 i = iper(j) 58 rmass(j) = rmass(i) 59c qres(j,:) = qres(i,:) 60 qres(j,1:idflow) = qres(i,1:idflow) 61 endif 62 enddo 63c 64c.... invert the diagonal mass matrix and find q 65c 66 rmass = one/rmass 67 68 do i=1, idflx 69 qres(:,i) = rmass*qres(:,i) 70 enddo 71 if (isurf .eq. 1) then 72 idflow=(nflow-1)*nsd 73c 74c.... calculation of the unit normal vector 75c 76 tmp = sqrt(qres(:,idflow+1)**2 77 & + qres(:,idflow+2)**2 78 & + qres(:,idflow+3)**2) 79 do i = 1, nshg 80 if (tmp(i) .lt. 0.0001) tmp(i) = 0.0001 81 end do 82 tmp = one/tmp 83 84 do i=1, nsd 85 qres(:,idflow+i) = tmp*qres(:,idflow+i) 86 enddo 87 endif 88 89 if(numpe > 1) then 90 call commu (qres, ilwork, idflx, 'out') 91 endif 92 93 if(iabc==1) then !are there any axisym bc's 94c 95c slave has masters value, for abc we need to rotate it 96c 97 do i=1,idflx/nsd 98 do j=1,2 99 istrt=j+(i-1)*(nflow-1) 100 uv(:,j)=qres(:,istrt) 101 enddo 102 call rotabc(uv, iBC, 'out') 103 do j=1,2 104 istrt=j+(i-1)*(nflow-1) 105 qres(:,istrt)=uv(:,j) 106 enddo 107 enddo 108 endif 109 110c 111c.... return 112c 113 return 114 end 115 116 117 subroutine qpbcSclr( rmass, qres, iBC, iper, ilwork ) 118c--------------------------------------------------------------------- 119c 120c This routine satisfies the periodic boundary conditions 121c on the diffusive flux residual and mass matrix 122c 123c input: 124c rmass (nshg) : mass matrix 125c qres (nshg, nsd) : diffusive flux vector 126c 127c output: modified qres and rmass 128c--------------------------------------------------------------------- 129 include "common.h" 130 131 dimension rmass(nshg), qres(nshg,nsd), 132 & iBC(nshg), iper(nshg) 133 134 if(iabc==1) !are there any axisym bc's 135 & call rotabc(qres, iBC, 'in ') 136 137c 138c.... compute qi for node A, i.e., qres <-- qres/rmass 139c 140 if (numpe > 1) then 141 call commu (qres , ilwork, nsd , 'in ') 142 call commu (rmass , ilwork, 1 , 'in ') 143 endif 144 145c 146c take care of periodic boundary conditions 147c 148 do j= 1,nshg 149 if (btest(iBC(j),10)) then 150 i = iper(j) 151 rmass(i) = rmass(i) + rmass(j) 152 qres(i,:) = qres(i,:) + qres(j,:) 153 endif 154 enddo 155 156 do j= 1,nshg 157 if (btest(iBC(j),10)) then 158 i = iper(j) 159 rmass(j) = rmass(i) 160 qres(j,:) = qres(i,:) 161 endif 162 enddo 163c 164c.... invert the diagonal mass matrix and find q 165c 166 rmass = one/rmass 167 168 do i=1, nsd 169 qres(:,i) = rmass*qres(:,i) 170 enddo 171 172 if(numpe > 1) then 173 call commu (qres, ilwork, nsd, 'out') 174 endif 175 176 if(iabc==1) !are there any axisym bc's 177 & call rotabc(qres, iBC, 'out') 178 179c 180c.... return 181c 182 return 183 end 184 185