FLOW VECTOR CANCELLING EACH OTHER FOR SOURCE/SINK TERMS
Dear,
I notice that for all the source and sink elements flow vector with opposite direction are cancelling each other. To avoid this I would like to output the absolute total flow, so instead of having in the code:
totflo(k,np) = flo(nnp+np)*dvecta(k,iconel) + totflo(k,np)
I was thinking to have something like:
otflo(k,np) = ABS(flo(nnp+np))*dvecta(k,iconel) + totflo(k,np)
but it doesn't work: 'The shapes of the array expressions do not conform. [TOTFLO]'. I have tried different things but nothing works.
Does anybody know how to solve this?
Thanks,
Pierre
6 replies

I forget to say that these are the total flow vector alonf the x,y and z axis outputted in the file flowvector.tec

Are you using the standard version of TOUGH3? I can not find the statement "totflo(k,np) = flo(nnp+np)*dvecta(k,iconel) + totflo(k,np)" in it.

Kenny Writing the fluxes are a special feature in TOUGHREACT/TReactMech. Yes, there are some issues on the fluxes from GENER which are not known in the fluxes at connection when the fluxes/area are calculated, particularly at boundaries. I'll take a look at it, since I created this! It gets messy because of the polygonal grids and the array variables to limit excess calculations and memory.

Hi Pierre,
I'm not exactly sure what you did here, but will take a closer look. The areaweighting is the dvecta factor, calculated in t2f in the following loop:
c..... Calculate all the geometric factors for the flux & velocity vector calculations
c... Fluxes are areaweighted, velocities not
!
!$omp parallel
!els012420 !$omp& shared (ncon,nex1,nex2,dvect,dvecta,x1,x2,x3)
!$omp& shared (ncon,nex1,nex2,dvect,dvecta,x1,x2,x3i,area)
!$omp& private (n,n1,n2,dxconn,dyconn,dzconn,vectdist)
!$omp DO SCHEDULE(static)
!
do n = 1, ncon
n1 = nex1(n)
n2 = nex2(n)
dxconn = x1(n1)x1(n2)
dyconn = x2(n1)x2(n2)
dzconn = x3(n1)x3(n2)
vectdist = 2.d0*sqrt(dxconn**2+dyconn**2+dzconn**2)
if(vectdist.eq.0.d0)vectdist = 1.d0
dvect(1,n) = dxconn/vectdist
dvecta(1,n) = dvect(1,n)/area(n)
dvect(2,n) = dyconn/vectdist
dvecta(2,n) = dvect(2,n)/area(n)
dvect(3,n) = dzconn/vectdist
dvecta(3,n) = dvect(3,n)/area(n)
enddo