Fix for node volume

An issue in the node volume calculation in cylindrical coordinates was
found. This was causing wrong conservation of current. Still to test
with ALPHIE_Grid case.

Still to check triangular element.
Still to theck 1D radial geometry
This commit is contained in:
Jorge Gonzalez 2022-12-20 15:51:43 +01:00
commit dd1fca3fee
2 changed files with 19 additions and 6 deletions

View file

@ -249,7 +249,7 @@ MODULE moduleMesh2DCyl
dr = self%r(2) - self%r(1)
dz = self%z(2) - self%z(1)
IF (dr /= 0.D0) THEN
r(2) = dr*DSQRT(rnd) + self%r(1)
r(2) = dr * DSQRT(rnd) + self%r(1)
r(1) = dz * (r(2) - self%r(1))/dr + self%z(1)
ELSE
@ -320,7 +320,18 @@ MODULE moduleMesh2DCyl
fPsi = self%fPsi(xi)
r = DOT_PRODUCT(fPsi,self%r)
self%volume = r*detJ
self%arNodes = fPsi*r*detJ
xi = (/-5.D-1, -5.D-1, 0.D0/)
r = DOT_PRODUCT(self%fPsi(xi),self%r)
self%arNodes(1) = fPsi(1)*r*detJ
xi = (/ 5.D-1, -5.D-1, 0.D0/)
r = DOT_PRODUCT(self%fPsi(xi),self%r)
self%arNodes(2) = fPsi(2)*r*detJ
xi = (/ 5.D-1, 5.D-1, 0.D0/)
r = DOT_PRODUCT(self%fPsi(xi),self%r)
self%arNodes(3) = fPsi(3)*r*detJ
xi = (/-5.D-1, 5.D-1, 0.D0/)
r = DOT_PRODUCT(self%fPsi(xi),self%r)
self%arNodes(4) = fPsi(4)*r*detJ
END SUBROUTINE areaQuad