Windkessel outlet boundary condition

Questions and remarks about code_saturne usage
Forum rules
Please read the forum usage recommendations before posting.
Saturne_User2021
Posts: 6
Joined: Sun Apr 18, 2021 4:44 pm

Windkessel outlet boundary condition

Post by Saturne_User2021 »

Hello there

is it possible to implement a three element Windkessel type outlet boundary condition with code_saturne?
To do so, I need to model a pressure-flow coupled boundary Condition. Here, a pressure outlet BC must be used, but its value depends on pressure and flow at previous time step (to approximate derivative terms).
How can I get access to these informations at the previous time step with a UDF?
Yvan Fournier
Posts: 4070
Joined: Mon Feb 20, 2012 3:25 pm

Re: Windkessel outlet boundary condition

Post by Yvan Fournier »

Hello,

I am not familiar with this type of boundary condition, but in any case, if you need values at the previous time steps, you could always create an array in a user module (or simply with a "save" attribute in cs_user_boundary_conditions.f90 or "static" in cs_user_boundary_conditions.c) to save values you need for the next time step.

In any case, when the boundary conditions UDF is called, we are preparing for a new time step, and pressure, velocity, ... values available are those of the previous time step, so you might actually not need to save anything, just access the data you need.

In case such an explicit boundary condition causes stability issues (I would not be surprised), some form of relaxation might be useful.

Best regards,

Yvan
Saturne_User2021
Posts: 6
Joined: Sun Apr 18, 2021 4:44 pm

Re: Windkessel outlet boundary condition

Post by Saturne_User2021 »

Hi Yvan

Thank you for your response.

If the UDF run is terminated before the current time step, there is definitely no need for data at the previous time step.
Anyways, I will test your recommendations, and I will keep you informed of progress. It may help other users of Code_Saturne in the medical fields.

Thank you for your assistance
Yvan Fournier
Posts: 4070
Joined: Mon Feb 20, 2012 3:25 pm

Re: Windkessel outlet boundary condition

Post by Yvan Fournier »

Hello,

Yes, the UDF run (for the boundary conditions UDF) is done before resolution of the new time step, so that should be good.

Best regards,

Yvan
Saturne_User2021
Posts: 6
Joined: Sun Apr 18, 2021 4:44 pm

Re: Windkessel outlet boundary condition

Post by Saturne_User2021 »

Hi Yvan,
I was digging in the doxygen document and I found the following modules to retrieve fields values at the preveious time step, e.g for pressure :

Code: Select all

! INITIALIZATION
call field_set_n_previous(ivarfl(ipr),1)
call field_get_val_prev_s(ivarfl(ipr), press_temp)
then

Code: Select all

pretem = press_temp(ifac)

Could you please confirm wether ("pretem") gives pressure at (-1) time step?

Thank you for all your help
Yvan Fournier
Posts: 4070
Joined: Mon Feb 20, 2012 3:25 pm

Re: Windkessel outlet boundary condition

Post by Yvan Fournier »

Hello,

Code: Select all

call field_set_n_previous(ivarfl(ipr),1)
Is not necessary here, as it is the default (and changing this for solved variables is risky, as it may not be necessary with a given time scheme, especially if you reduce the value).
call field_get_val_prev_s(ivarfl(ipr), press_temp)
Will get the value at time step n-1, but if you call it while preparing boundary conditions for time step n+1, it is probably not the value you want. You may prefer
call field_get_val_s(ivarfl(ipr), press_temp)
And finally, the provides a cell-based field, not a face-based one. You can get the cell if matching a given boundary face using "ifabor(ifac)" in Fortran, or "cs_glob_mesh->b_face_cells[face_id]" in C.

If you really need to face values instead of matching cell values (to avoid interpolation issues), you can use the non-reconstructed boundary value using a First-order Taylor expansion by also calling a gradient operator, but for starts, I would try with the simple cell value, which should be very close.

Best regards,

Yvan
Saturne_User2021
Posts: 6
Joined: Sun Apr 18, 2021 4:44 pm

Re: Windkessel outlet boundary condition

Post by Saturne_User2021 »

Thank you Yvan for detailed explanation.
Since I am dealing with derivative terms (d/dt) for both pressure and flow rate, I need at least two values at current and previous time steps to set the BC. That is the reason I am trying to get access to previous values ( t =0 and -1 or -2 here).

It seems to be working fine for the pressure, however, I get an out of bound error using te following lines to get the mass flux at previous time :

Code: Select all

call field_get_key_int(ivarfl(ivar), kimasf, iflmas)
call field_get_key_int(ivarfl(ivar), kbmasf, iflmab)

call field_set_n_previous(iflmas,1)
call field_set_n_previous(iflmab,1)

call field_get_val_prev_s(iflmas, imasfl_prev)
call field_get_val_prev_s(iflmab, bmasfl_prev)
The error message I get :
SIGSEGV signal (forbidden memory area access) intercepted!

Call stack:
1: 0x7feec1bc1760 <cs_restart_write_section+0x1d0> (libsaturne.so.5)
2: 0x7feec1bce7ca <cs_restart_write_field_vals+0x7a> (libsaturne.so.5)
3: 0x7feec1bceb3f <cs_restart_write_linked_fields+0x33f> (libsaturne.so.5)
4: 0x7feec21fcea7 <__cs_c_bindings_MOD_restart_write_linked_fields+0xfc> (libsaturne.so.5)
5: 0x7feec1c123d6 <ecrava_+0x14d8> (libsaturne.so.5)
6: 0x7feec1af9004 <caltri_+0x1fb5> (libsaturne.so.5)
7: 0x7feec32854ff <cs_run+0x57f> (libcs_solver.so.5)
8: 0x7feec3284e35 <main+0x125> (libcs_solver.so.5)
9: 0x7feec0ce2b97 <__libc_start_main+0xe7> (libc.so.6)
10: 0x557b51e4f02a <_start+0x2a> (cs_solver)
End of stack
Any suggestion here?
Yvan Fournier
Posts: 4070
Joined: Mon Feb 20, 2012 3:25 pm

Re: Windkessel outlet boundary condition

Post by Yvan Fournier »

Hello,

Yes, by default, the mass flux is only stored for 1 time step. Forcing it to 2 time steps should probably work, in which case you can use :

Code: Select all

call field_set_n_previous(iflmas,1)
call field_set_n_previous(iflmab,1)
But this must be called before the time loop, for example in cs_user_parameters.f90, in the "usipsu" function.
And it is safer not to call this again in the time loop (i.e. for the boundary conditions), but simply to access the values..

Best regards,

Yvan
Saturne_User2021
Posts: 6
Joined: Sun Apr 18, 2021 4:44 pm

Re: Windkessel outlet boundary condition

Post by Saturne_User2021 »

Hello Yvan

As suggested in your posting, I added the following lines to cs_user_parameters.f90.

Code: Select all

subroutine usipsu &
 ( nmodpp )

!===============================================================================
! Module files
!===============================================================================

use paramx
use cstnum
use dimens
use numvar
use optcal
use cstphy
use entsor
use parall
use period
use ihmpre
use albase
use ppppar
use ppthch
use ppincl
use coincl
use cpincl
use field
use cavitation
use post
use rotation
use cs_c_bindings

!===============================================================================

implicit none

! Arguments
integer nmodpp

integer iflmas,iflmab
!===============================================================================

!>  This subroutine allows setting parameters
!>  which do not already appear in the other subroutines of this file.
!>
!>  It is possible to add or remove parameters.
!>  The number of physical properties and variables is known here.

!===============================================================================

!----
! Formats
!----

call field_set_n_previous(iflmas,1)
call field_set_n_previous(iflmab,1)

return
end subroutine usipsu
But the simulation crashed with the following message :

/home/A41771/salome/V8_5_0/modules/src/PRECFDSTUDY/src/base/cs_field.c:2279: Fatal error.

Field with id -586798064 is not defined.


Call stack:
1: 0x7f70db526025 <cs_field_by_id+0x55> (libsaturne.so.5)
2: 0x7f70db526148 <cs_f_field_set_n_previous+0x8> (libsaturne.so.5)
3: 0x7f70dbbaf93b <__field_MOD_field_set_n_previous+0xd> (libsaturne.so.5)
4: 0x55680df71c2f <usipsu_+0x15> (cs_solver)
5: 0x7f70db5dd616 <iniusi_+0x3ac> (libsaturne.so.5)
6: 0x7f70db5dce8a <initi1_+0x63> (libsaturne.so.5)
7: 0x7f70dcc3f44e <cs_run+0x4ce> (libcs_solver.so.5)
8: 0x7f70dcc3ee35 <main+0x125> (libcs_solver.so.5)
9: 0x7f70daa7bb97 <__libc_start_main+0xe7> (libc.so.6)
10: 0x55680df71b3a <_start+0x2a> (cs_solver)
End of stack
It seems to be tough to force the mass flux valued to -1. If so, do you know any similar method to obtain the flow-rate at previous time step? Maybe with normal vectors and velocity components ?
Thanks a lot for your help
Yvan Fournier
Posts: 4070
Joined: Mon Feb 20, 2012 3:25 pm

Re: Windkessel outlet boundary condition

Post by Yvan Fournier »

Hello,

Sorry, I forgot iflmas/iflmab may not be initialed at that point.

If you move the initialization to usipes (in the same file), which is called a bit later, things should be better.

Best regards,

Yvan
Post Reply