0

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

null
    • Pierre_Jeanne
    • 11 mths ago
    • Reported - view

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

    • kenny
    • 11 mths ago
    • Reported - view

    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. 

    • TOUGHREACT Developer
    • Eric_Sonnenthal
    • 11 mths ago
    • Reported - view

    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.

      • Pierre_Jeanne
      • 11 mths ago
      • Reported - view

       

      Hello Kenny and Eric,

      thank you for taking a look at that. I added a few lines to output  the flow  along all the connections associated with an element. I am happy with that, but it will be even better if I could also normalized them to their respective area.. I didn't find away to do that yet...

              do m = 1, nconel
                iconel = icelem(m,i)          
                nnp = (iconel-1)*nph
                np = 1
                cflo(m,np) = ABS(FLO(nnp+np))     ---------- I added this ---
                   do k = 1, 3
      c.... Fluid flux normalized to area of connection
                    totflo(k,np) = flo(nnp+np)*dvecta(k,iconel) + totflo(k,np)
                  enddo   
              enddo

       

      Thank you!

    • TOUGHREACT Developer
    • Eric_Sonnenthal
    • 11 mths ago
    • Reported - view

    Hi Pierre,

     

    I'm not exactly sure what you did here, but will take a closer look. The area-weighting 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 area-weighted, 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

      • Pierre_Jeanne
      • 11 mths ago
      • Reported - view

       

      Yes, you were right what I did does not work...

      the issue that I have is that for example if I have a horizontal 2D model and i simulate injection in  one element present in a fracture oriented along the x-axis in an homogeneous reservoir then the flow in the y axis equal 0 in all the grid element composing the fracture. The flows have opposite direction and must cancel each other.

Content aside

  • 11 mths agoLast active
  • 6Replies
  • 35Views
  • 3 Following