Incorrect value for convective outlet boundary condition coefficients in condli.f90
Posted: Sun Aug 30, 2020 11:53 am
Hello,
As it is mentioned in the theory PDF file of Code_Saturne 6.0.4, for the convective outlet boundary condition we have:
As I noticed the A corresponds to coefa and B corresponds to coefb in the condli.f90 file.
But in that file for the set_convective_outlet_vector subroutine, the value of B is:
Which means:
This mistake has also been taken place in the set_convective_outlet_tensor and set_convective_outlet_vector_aniso subroutines.
========================================
Also the input values for this type of B.C was confusing for me and some users, so I checked the code and found out that for this type of B.C we have:
rcodcl(ifac,ivar,1) = The value of ivar at previous time step
rcodcl(ifac,ivar,2) = CFL number based on the celerity value.
Hope that this will be helpful and will be corrected in the next versions.
Regards,
Mohammad
As it is mentioned in the theory PDF file of Code_Saturne 6.0.4, for the convective outlet boundary condition we have:
I checked these values by discretizing the equation and they are 100% correct.A = (1/(1+CFL))Y
B = (CFL/(1+CFL))
As I noticed the A corresponds to coefa and B corresponds to coefb in the condli.f90 file.
But in that file for the set_convective_outlet_vector subroutine, the value of B is:
Code: Select all
coefb(isou,jsou) = cflv(isou)*(1.d0+cflv(isou))
Now you can see that the "*" must be replaced by "/" in the code which leads to:B = (CFL*(1+CFL))
Code: Select all
coefb(isou,jsou) = cflv(isou)/(1.d0+cflv(isou))
========================================
Also the input values for this type of B.C was confusing for me and some users, so I checked the code and found out that for this type of B.C we have:
rcodcl(ifac,ivar,1) = The value of ivar at previous time step
rcodcl(ifac,ivar,2) = CFL number based on the celerity value.
Hope that this will be helpful and will be corrected in the next versions.
Regards,
Mohammad