BDF coefficient vector length in FluidDyn app

[pablobecker]

Hi, I was getting random with an element created with the symbolic generation, so the header file is similar to others of the fluid dynamic app.

Finally I was able to track the problem to the BDF coefficients vector. The problem is:

• As other python fluid solvers, we are using the ComputeBDFCofficients to set this vector. (e.g. https://github.com/KratosMultiphysics/Kratos/blob/master/applications/FluidDynamicsApplication/python_scripts/navier_stokes_compressible_solver.py#L132 )
• For stability reasons, we are using BDF1, time order =1
• This means that the BDF vector is of size 2: https://github.com/KratosMultiphysics/Kratos/blob/master/kratos/processes/compute_bdfcoefficients_process.h#L149 *However the fluid element is reading the BDF vector as it had size=3 without even checking https://github.com/KratosMultiphysics/Kratos/blob/master/applications/FluidDynamicsApplication/custom_elements/compressible_navier_stokes.h#L402

Hence the random results.

I was about to create a PR to modify the affected header elements so that the size is first checked, but decided to create this issue so that we discuss what is the best way to fix it. Let me know what you think @RiccardoRossi @jcotela @rubenzorrilla

[rubenzorrilla]

So far, all the symbolic elements assume that the BDF2 is used and thus that the BDF coefficients vector has size 3. I wouldn't be in favor of introducing an if statement to handle this.

In any case, the problem comes from the fact that the symbolic generation is quite rigid in this aspect, and requires to known the acceleration formula (what is to say the time scheme) a priori to check which terms of the acceleration need to be differentiated.

[RiccardoRossi]

well we will need to fix this somehow...it should not be difficult to do it in the template

On Thu, Jul 4, 2019, 1:59 PM Rubén Zorrilla notifications@github.com wrote:

So far, all the symbolic elements assume that the BDF2 is used and thus that the BDF coefficients vector has size 3. I wouldn't be in favor of introducing an if statement to handle this.

In any case, the problem comes from the fact that the symbolic generation is quite rigid in this aspect, and requires to known the acceleration formula (what is to say the time scheme) a priori to check which terms of the acceleration need to be differentiated.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/KratosMultiphysics/Kratos/issues/5175?email_source=notifications&email_token=AB5PWEP3BTVCLYSWRNEZZHTP5XX4DA5CNFSM4H5ZZHXKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODZHLFZI#issuecomment-508474085, or mute the thread https://github.com/notifications/unsubscribe-auth/AB5PWENS7AJCTXZNAQY2CS3P5XX4DANCNFSM4H5ZZHXA .

[pablobecker]

the solver settings allow to set it to bdf1, so it is accepting it

[pablobecker]

moreover, IMO it is important to allow to use BDF1. with variable timesteps it helps a lot in stability

[RiccardoRossi]

anyhow, we are also accessing the velocities and pressures two steps back...

one option woukd be to keep a buffer size of 3 and only use the first two components (third component set to zero)

On Thu, Jul 4, 2019, 2:04 PM Pablo notifications@github.com wrote:

the solver settings allow to set it to bdf1, so it is accepting it

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/KratosMultiphysics/Kratos/issues/5175?email_source=notifications&email_token=AB5PWENESI7VREFVDOLDIC3P5XYNLA5CNFSM4H5ZZHXKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODZHLRAI#issuecomment-508475521, or mute the thread https://github.com/notifications/unsubscribe-auth/AB5PWEIDALCKECO42TL34STP5XYNLANCNFSM4H5ZZHXA .

[jcotela]

If I understand this correctly, the symbolic element, once compiled, assumes a given BDF order, so the size of the vector is so far "fixed" the moment the symbolic code is generated.

Ideally, we would always check for the size of the BDF vector during the corresponding Element::Check call. I remember having tried this a long time ago but having had to roll it back due to the Check method being called before the process that initializes the BDF. I do not know if the changes in the Solver in the meantime (this was before the AnalysisStage existed) have helped in this situation, but maybe this is something worth keeping in mind.

the solver settings allow to set it to bdf1, so it is accepting it

The element you pointed us to is pretty much assuming BDF2... If it is accepting BDF1 I would be inclined to call that a bug. The situation may be different for your own element, of course, but then you have the choice to support it on the level of the specific element.

[rubenzorrilla]

If I understand this correctly, the symbolic element, once compiled, assumes a given BDF order, so the size of the vector is so far "fixed" the moment the symbolic code is generated.

Ideally, we would always check for the size of the BDF vector during the corresponding Element::Check call. I remember having tried this a long time ago but having had to roll it back due to the Check method being called before the process that initializes the BDF. I do not know if the changes in the Solver in the meantime (this was before the AnalysisStage existed) have helped in this situation, but maybe this is something worth keeping in mind.

the solver settings allow to set it to bdf1, so it is accepting it

The element you pointed us to is pretty much assuming BDF2... If it is accepting BDF1 I would be inclined to call that a bug. The situation may be different for your own element, of course, but then you have the choice to support it on the level of the specific element.

Exactly, accepting in this case is actually a bug (possibly my bad)

[pablobecker]

But as I pointed in the first link ( https://github.com/KratosMultiphysics/Kratos/blob/master/applications/FluidDynamicsApplication/python_scripts/navier_stokes_compressible_solver.py#L132 ) , it has settings to change it.

Although with a different script, I am also using the stokes element of this App. It works OK with BDF1 if BDF[2] =0 . We are accessing unnecesarily older velocities and pressure, but nevertheless I'd say it is supported.

[rubenzorrilla]

Which element are you are aiming to use? The compressible navier-stokes element or the stokes one? This is important since the error comes from the element Automatic Differentiation (AD). As you point, the solver and the scheme "supports" it (even though the second order is hardcoded by default).

[pablobecker]

Stokes two fluid

[RiccardoRossi]

so let s resize the bdfvector to 3 and use a buffer skze of 3. it is suboptimal but jt eill work

On Thu, Jul 4, 2019, 2:25 PM Pablo notifications@github.com wrote:

Stokes two fluid

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/KratosMultiphysics/Kratos/issues/5175?email_source=notifications&email_token=AB5PWELOBKG26KF7X7PAKATP5X24XA5CNFSM4H5ZZHXKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODZHNCVA#issuecomment-508481876, or mute the thread https://github.com/notifications/unsubscribe-auth/AB5PWEISIW4J3SV5W2X4ZCLP5X24XANCNFSM4H5ZZHXA .

[pablobecker]

That would work. Assuming noone is using the size of this vector to determine time order

[pablobecker]

ping @jcotela @rubenzorrilla

[jcotela]

@pablobecker I see two options, either you implement your element to detect how many steps you want to use (by the size of the BDF vector or otherwise) or you make it so that the BDF vector always has size 3 (and add a 0 if needed for first order), for example by implementing your own BDF process.

Having said that, I think that the approach that you and @RiccardoRossi suggest is good only as a short-term solution and I would implement it only as an ad-hoc fix for your specific solver. I do not believe that keeping an extra time step (meaning 50% more nodal data memory in this case) to support multiple time orders is something that should be encouraged by giving it standard support across the application.

[pablobecker]

@jcotela . Yes, I see your point. As a workaround, we are currently using the modified BDF process with size3. Eventually we'll create first order elements in our apps and also for the stokes two phase of the fluiddynapp.

In any case, some python solvers of the fluiddyn app are wrong bc they accept BDF1, while the element does not. What do you want me to do, do I push the temporary workaround of the BDF process or you fix the python solver, allowing only BDF2 ?

[RiccardoRossi]

a nicer fix is to adapt this code:

https://github.com/KratosMultiphysics/Kratos/blob/master/applications/FluidDynamicsApplication/custom_elements/stokes_3D_twofluid.h#L162

here we can say that if the time order is one than we set to zero the variables IN THE LOCAL MATRIX.

On Tue, Jul 9, 2019 at 12:44 PM Pablo notifications@github.com wrote:

@jcotela https://github.com/jcotela . Yes, I see your point. As a workaround, we are currently using the modified BDF process with size3. Eventually we'll create first order elements in our apps and also for the stokes two phase of the fluiddynapp.

In any case, some python solvers of the fluiddyn app are wrong bc they accept BDF1, while the element does not. What do you want me to do, do I push the temporary workaround of the BDF process or you fix the python solver, allowing only BDF2 ?

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/KratosMultiphysics/Kratos/issues/5175?email_source=notifications&email_token=AB5PWENOYAXFY2BNNMF5NYDP6RT2PA5CNFSM4H5ZZHXKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODZP4CZY#issuecomment-509591911, or mute the thread https://github.com/notifications/unsubscribe-auth/AB5PWELNCOZ6QDEEVEROKULP6RT2PANCNFSM4H5ZZHXA .

--

Riccardo Rossi

PhD, Civil Engineer

member of the Kratos Team: www.cimne.com/kratos

Associate Professor at Universitat Politècnica de Catalunya, BarcelonaTech (UPC)

Full Research Professor at International Center for Numerical Methods in Engineering (CIMNE)

C/ Gran Capità, s/n, Campus Nord UPC, Building C1, First Floor

08034 – Barcelona – Spain – www.cimne.com -

T.(+34) 93 401 56 96 skype: rougered4

http://www.cimne.com/

https://www.facebook.com/cimne http://blog.cimne.com/ http://vimeo.com/cimne http://www.youtube.com/user/CIMNEvideos http://www.linkedin.com/company/cimne https://twitter.com/cimne

Les dades personals contingudes en aquest missatge són tractades amb la finalitat de mantenir el contacte professional entre CIMNE i voste. Podra exercir els drets d'accés, rectificació, cancel·lació i oposició, dirigint-se a cimne@cimne.upc.edu. La utilització de la seva adreça de correu electronic per part de CIMNE queda subjecte a les disposicions de la Llei 34/2002, de Serveis de la Societat de la Informació i el Comerç Electronic.

Imprimiu aquest missatge, només si és estrictament necessari. http://www.cimne.com/

[pablobecker]

to do that the element should know the time order we are using. Is any of the other apps saving the time order in the process info or something else? something like TIME_ORDER

[RiccardoRossi]

you could get it from the lenght of the BDF_VECTOR as you mentioned.

essentially what i am proposing is to do something like

  bdf_coeff0 = BDFVector[0]
  bdf_coeff1 = BDFVEctor[1]
  if(BDFVector.size() == 3)
       bdf_coeff2 = BDFVEctor[2]
  else
       bdf_coeff2 = 0.0

same when we access to the velocity two steps back

[philbucher]

Is this solved?

[rubenzorrilla]

I don't think so... I'd like to keep it as a reminder. I'd add it to the application project too.