Axail fan modelling usin cs_fan.h

Hello,

I’m trying to modelize an axial fan using the cs_fan.h.

Looking at the doxygen documentation in the programmer’s documentation, things seem to have change since http://code-saturne.org/forum/viewtopic.php?f=2&t=1045&hilit=axial+fan

As far as I’ve understood things :
In cs_user_initialization.f90, I’ve to define and initialze the fan.
In cs_user_source_terms, I’ve to tell to CS to compute the added source terms with tsvvtl.

For the moment, i’ve added those lines ine cs_user_initialization-fan.f90 in the SRC directory of the case, after the [init] :

 ! Fan initialization
!
!



call cs_fan_define(1, 3, (/ 0., 0., 0. /), (/ 0., 0., 0.2 /), .707d0, .65d0, .1d0, (/ 300., -14.67, 0. /), 1.d0)

! 1 : constant_f
! 3 : fan dimension 3 : 3D
! cu : intersection coords. of axis and inlet face
! cd : intersection coords. of axis and outlet face
! .707d0 : fan radius
! .65d0 : blades radius
! .1d0 : hub radius
! fcar : curve coefficients, dgre 0,1 and 2 of the characteristic curve
! 1.d0 : fan axial torque



!call inivtl()
! call numvtl(fannbr)
! print *, fannbr

When trying to run a calculation, I get the followings in the compile.log :

cygwin warning:
  MS-DOS style path detected: C:\tmp\TestCodeSaturne\Fan\fan1\RESU\20160222-1829\src_saturne\cs_user_initialization-fan.f90
  Preferred POSIX equivalent is: /cygdrive/c/tmp/TestCodeSaturne/Fan/fan1/RESU/20160222-1829/src_saturne/cs_user_initialization-fan.f90
  CYGWIN environment variable option "nodosfilewarning" turns off this warning.
  Consult the user's guide for more details about POSIX paths:
    http://cygwin.com/cygwin-ug-net/using.html#using-pathnames
C:\tmp\TestCodeSaturne\Fan\fan1\RESU\20160222-1829\src_saturne\cs_user_initialization-fan.f90:48.5:

   dt     )
     1
Warning: Unused dummy argument 'dt' at (1)
C:\tmp\TestCodeSaturne\Fan\fan1\RESU\20160222-1829\src_saturne\cs_user_initialization-fan.f90:92.28:

integer          iel, iutile
                            1
Warning: Unused variable 'iutile' declared at (1)
C:\tmp\TestCodeSaturne\Fan\fan1\RESU\20160222-1829\src_saturne\cs_user_initialization-fan.f90:47.17:

 ( nvar   , nscal  ,                                              &
                 1
Warning: Unused dummy argument 'nscal' at (1)
C:\tmp\TestCodeSaturne\Fan\fan1\RESU\20160222-1829\src_saturne\cs_user_initialization-fan.f90:47.7:

 ( nvar   , nscal  ,                                              &
       1
Warning: Unused dummy argument 'nvar' at (1)
cygwin warning:
  MS-DOS style path detected: C:\tmp\TestCodeSaturne\Fan\fan1\RESU\20160222-1829\cs_solver.exe
  Preferred POSIX equivalent is: /cygdrive/c/tmp/TestCodeSaturne/Fan/fan1/RESU/20160222-1829/cs_solver.exe
  CYGWIN environment variable option "nodosfilewarning" turns off this warning.
  Consult the user's guide for more details about POSIX paths:
    http://cygwin.com/cygwin-ug-net/using.html#using-pathnames
cs_user_initialization-fan.o:cs_user_initialization-fan.f90:(.text+0x15f): undefined reference to `cs_fan_define_'
collect2: error: ld returned 1 exit status
"C:\Program Files (x86)\Code_Saturne\4.0\bin\x86_64-w64-mingw32-gfortran" -I C:\tmp\TestCodeSaturne\Fan\fan1\RESU\20160222-1829\src_saturne -I "C:\Program Files (x86)\Code_Saturne\4.0\include\code_saturne" -x f95-cpp-input -Wall -pedantic-errors -std=f2003 -O -c C:\tmp\TestCodeSaturne\Fan\fan1\RESU\20160222-1829\src_saturne\cs_user_initialization-fan.f90
"C:\Program Files (x86)\Code_Saturne\4.0\bin\x86_64-w64-mingw32-g++" -o C:\tmp\TestCodeSaturne\Fan\fan1\RESU\20160222-1829\cs_solver.exe cs_user_initialization-fan.o "-LC:\Program Files (x86)\Code_Saturne\4.0\lib" -L/home/MFEE/opt/x86_64-w64-mingw32/lib -L/home/MFEE/opt/x86_64-w64-mingw32/lib -L/home/MFEE/opt/x86_64-w64-mingw32/lib -L/home/MFEE/opt/x86_64-w64-mingw32/lib -L/home/MFEE/opt/MPI/lib -L/home/MFEE/opt/x86_64-w64-mingw32/lib -O -Wl,--allow-multiple-definition -lsaturne -lple -lcgns -lmedC -lhdf5 -lscotch -lmsmpi -lxml2 -lz -lz -L/usr/lib/gcc/x86_64-w64-mingw32/4.8.2 -L/usr/lib/gcc/x86_64-w64-mingw32/4.8.2/../../../../x86_64-w64-mingw32/lib/../lib -L/usr/x86_64-w64-mingw32/sys-root/mingw/lib/../lib -L/usr/lib/gcc/x86_64-w64-mingw32/4.8.2/../../../../x86_64-w64-mingw32/lib -L/usr/x86_64-w64-mingw32/sys-root/mingw/lib -lz -lgfortran -lmingw32 -lmoldname -lmingwex -lmsvcrt -lquadmath -ladvapi32 -lshell32 -luser32 -lkernel32 -lm

Does anybody know what could be wrong?

As you may have found, I’m running on Windows 7 64 bits
Code_saturne version is 4.0.2

Thanks a lot for your help…

JB

Hello,

Actually, you need to call fans from C code (cs_user_parameters.c instead of cs_user_parameters.f90).

Unfortunately, I don’t believe fans work well in 4.0: they were an update of an old, unfinished development, which was not much used, and not part of the validation suite. The fan functionality has seen new use and testing very recently, and is in a mostly working state starting with version 4.2, with some fixes in the branches but not yet in releases (so you would need to add cs_fan.c from the trunk: http://code-saturne.org/viewvc/saturne/trunk/src/base/cs_fan.c?revision=8520&view=markup, as the fix has not been merged to branch 4.2 yet; i’ll try to remember to merge it tomorrow).

In any case, I would not recommend using this with any version older than 4.2…

Note that the theory documentation for fan is not available yet: I have an old, partial documentation from an older code in which that functionnality existed, but we need to translate it from French ton English and from a pdf to LaTeX or Doxygen…

Regards,

Yvan

thanks a lot for your answer, Yvan.

As I run code_saturne on Windows, is there any recompilation to do, or I just have to replace the cs_fan.c file in the installation directory by the one found in the trunk?

By the way, do you know if it is possible to have two Code_saturne installation in parallel on windows? A 4.0.3 and a 4.2? Or is it better to have only one? In this case, I’ll have to update my install.

Jb

Hello,

Adding the cs_fan.c from the head of the 4.2 branch to your SRC directory (assuming you switch to version 4.2) should work.

Installing 2 Code_Saturne versions in parallel works well on Linux, and should probably work on Windows (though I need to check with Cyril on that), but I am not sure it has been tested.

Regards,

Yvan

Hello,

I’ve tried to copy cs_fan.h, cs_fan.c and cs_user_parameters-fan.c from a CS 4.2.1 installation on Linux in the SRC directory of a case build with CS 4.0.3 on Windows.
The compilation fails. This way to do doesn’t seems to be very clean.

So, in order to overcome this issue, I’ve try to modelize a simple case of an internal fan with CS 4.2.1 on Linux.

The mesh was built in Salome.
The domaine is a box of 0.5 m side length
The mesh is an hexaedral one, with cell size of 0.02m in every direction.
External faces of the domaine have been grouped in a “Wall” face group.

I’ve then set my case in CS, using stationnary simulation and air as fluid.
I’ve copy cs_user_parameter-fan.c in the SRC directory, and changing the fan setting as follow :

cs_user_parameters(void)
{
  /*
   * We define a fan, which will be handled as an automatic
   * explicit source term. See cs_fan_define() prototype
   * or Doxygen documentation for details.
   */

  /*! [fan_user_1] */
  {
    const cs_real_t  inlet_axis_coords[3] = {0.22, 0.25, 0.25};
    const cs_real_t  outlet_axis_coords[3] = {0.28, 0.25, 0.25};
    const cs_real_t  fan_radius = 0.05;
    const cs_real_t  blades_radius = 0.02;
    const cs_real_t  hub_radius = 0.02;
    const cs_real_t  pressure_curve_coeffs[3] = {0.6, -0.1, -0.05};
    const cs_real_t  axial_torque = 0.001;

    cs_fan_define(3, /* fan (mesh) dimension */
                  inlet_axis_coords,
                  outlet_axis_coords,
                  fan_radius,
                  blades_radius,
                  hub_radius,
                  pressure_curve_coeffs,
                  axial_torque);
  }
  /*! [fan_user_1] */
}

Then, I ran the simulation.
There is no compiling error and simulation run.

But when checking the listing file, it appears that there is no fluid motion in the domain.
Velocity stay at 0! :astonished:

What did I do wrong in my set up?

Please find enclosed the listing file.

Thanks a lot for your help.

JB
listing.txt (73.6 KB)

With a reread of the code of cs_user_parameter-fan.c posted above, I’ve seen the following mistakes :

Blade radius = 0.02 and hub radius = 0.02 too.
So, there is no flow for a fan without blade…
Damnely true…

Moreover, the pressure curve coeffs aren’t very phisical for my case.
So, I’ve changed them to something found for an electronic cooling fan of 92mm diameter.

The cs_user_parameter-fan.c is now :

cs_user_parameters(void)
{
  /*
   * We define a fan, which will be handled as an automatic
   * explicit source term. See cs_fan_define() prototype
   * or Doxygen documentation for details.
   */

  /*! [fan_user_1] */
  {
    const cs_real_t  inlet_axis_coords[3] = {0.22, 0.25, 0.25};
    const cs_real_t  outlet_axis_coords[3] = {0.28, 0.25, 0.25};
    const cs_real_t  fan_radius = 0.05;
    const cs_real_t  blades_radius = 0.048;
    const cs_real_t  hub_radius = 0.02;
    const cs_real_t  pressure_curve_coeffs[3] = {40.88, -1474.7, 4587.1};
    const cs_real_t  axial_torque = 0.001;

    cs_fan_define(3, /* fan (mesh) dimension */
                  inlet_axis_coords,
                  outlet_axis_coords,
                  fan_radius,
                  blades_radius,
                  hub_radius,
                  pressure_curve_coeffs,
                  axial_torque);
  }
  /*! [fan_user_1] */
}

The subroutines compile and the calculation run until the end, but the listing file is very strange.
From iteration 3, it’s look like diverging :

INSTANT    0.300000000E-01   TIME STEP NUMBER               3
 ============================================================= 




   ** BOUNDARY MASS FLOW INFORMATION
      ------------------------------

---------------------------------------------------------------
Boundary type          Code    Nb faces           Mass flow
---------------------------------------------------------------
Inlet                     2           0         0.000000000E+00
Smooth wall               5        3750         0.000000000E+00
Rough wall                6           0         0.000000000E+00
Symmetry                  4           0         0.000000000E+00
Free outlet               3           0         0.000000000E+00
Free inlet               13           0         0.000000000E+00
Convective inlet         15           0         0.000000000E+00
Free surface             14           0         0.000000000E+00
Undefined                 1           0         0.000000000E+00
---------------------------------------------------------------


   ** BOUNDARY CONDITIONS FOR SMOOTH WALLS
   ---------------------------------------
------------------------------------------------------------
                                         Minimum     Maximum
------------------------------------------------------------
   Rel velocity at the wall uiptn :  0.00000E+00 0.00000E+00
   Friction velocity        uet   :  0.00000E+00 0.00000E+00
   Friction velocity        uk    :  0.10000E+13-0.10000E+13
   Dimensionless distance   yplus :  0.10000E+13-0.10000E+13
   ------------------------------------------------------
   Nb of reversal of the velocity at the wall   :          0
   Nb of faces within the viscous sub-layer     :       3750
   Total number of wall faces                   :       3750
------------------------------------------------------------


@                                                            
@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
@                                                            
@ @@ WARNING: MESH TOO REFINED AT THE WALL                   
@    ========                                                
@    The mesh seems to be too refined at the wall to use     
@      a high-Reynolds turbulence model.                     
@                                                            
@    The last time step at which too small values for the    
@      dimensionless distance to the wall (yplus) have been  
@      observed is the time step          3
@                                                            
@    The minimum value for yplus must be greater than the    
@      limit value YPLULI =    0.23810E+01
@                                                            
@    Have a look at the distribution of yplus at the wall    
@      (with EnSight or ParaView for example) to conclude on 
@      the way the results quality might be affected.        
@                                                            
@    This warning is only printed at the first two           
@      occurences of the problem and at the last time step   
@      of the calculation. The vanishing of the message does 
@      not necessarily mean the vanishing of the problem.    
@                                                            
@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
@                                                            


   ** INFORMATION ON CONVERGENCE
      --------------------------
-----------------------------------------------------------------------------
   Variable    Rhs norm      N_iter  Norm. residual      Drift  Time residual
-----------------------------------------------------------------------------
c  Velocity             NaN       1   0.00000E+00           NaN         NaN 
c  Velocity[X]                                              NaN 
c  Velocity[Y]                                              NaN 
c  Velocity[Z]                                              NaN 
c  Pressure             NaN       0   0.00000E+00           NaN 0.00000E+00 
c  k                    NaN       1   0.00000E+00           NaN         NaN 
c  epsilon              NaN       1   0.00000E+00           NaN         NaN 
-----------------------------------------------------------------------------
  
  

  ** "Clippings" des champs calculés
     -------------------------------

   champ                min. init.      max. init.   lim. à min.   lim. à max.
-  ---------------  --------------  --------------  ------------  ------------
a  k              
a  epsilon        
a  Local Time Step


  ** Champs calculés sur cellules
     ----------------------------

   champ                   minimum         maximum   moy. ensemble   moy. spatiale
-  ---------------  --------------  --------------  --------------  --------------
v  Velocity[X]                 inf            -inf            -nan            -nan
v  Velocity[Y]                 inf            -inf            -nan            -nan
v  Velocity[Z]                 inf            -inf            -nan            -nan
v  ǁVelocityǁ                  inf            -inf            -nan            -nan
v  Pressure                      0               0               0               0
v  k                           inf            -inf            -nan            -nan
v  epsilon                     inf            -inf            -nan            -nan
v  TurbVisc                    inf            -inf            -nan            -nan
v  CourantNb                   inf            -inf            -nan            -nan
v  FourierNb                   inf            -inf            -nan            -nan
v  total_pressure       1.0132e+05      1.0132e+05      1.0132e+05      1.0132e+05
v  Local Time Step       0.0098945         0.01331        0.013308        0.013308

  ** Champs calculés sur faces de bord
     ---------------------------------

   champ         minimum         maximum   moy. ensemble   moy. spatiale
-  -----  --------------  --------------  --------------  --------------
v  Yplus             inf            -inf            -nan            -nan


  ** Valeurs calculées sur cellules
     ------------------------------

   criterion               minimum         maximum   moy. ensemble   moy. spatiale
   ---------------  --------------  --------------  --------------  --------------
   Courant/Fourier             inf            -inf            -nan            -nan


Fans
----
    id      surface       volume         flow       deltaP
  ----  -----------  -----------  -----------  -----------
     0   4.0800e-02   8.0000e-05         -nan        -nan

Does anybody know what is wrong?

Please find enclosed the full listing

JB
listing.txt (2 MB)

Hello,

The code actually seems to diverge from iteration 1…

Could you try running on a single rank (i.e. non-parallel) for an iteration or 2 ? I don’t expect this will change anything, but it is a quick check, so I recommend doing it “just in case”.

Also, do you have access to or can you install a debug build ? A debug build is slower, but should catch the divergence (or division by zero, …) much closer to the cause of the error, making it much easier to analyze and fix.

Finally, increasing verbosity to 2 for all variables (search for “iwarni” in the user subroutine examples) will also help pinpoint where the problem originates.

Best regards,

Yvan

Hello,

Unfortunately, I don’t have access to a debug build.
The solution would be to install another instance of CS with debug option on my Linux install.
Why not, but as it’s quite complicated for me to properly install CS, I would like to explore other options before.

I’ve tried the same case using serial solver.
This doesn’t seem to change anything.

I’ve also put IWARNI variables to 2 and relaunch the simulation.
Here is the listing file of the run.

Something is strange to me : why there is only 1 iteration for velocity, k and epsilon?
Also, the maximum number of iteration seems to have been removed from the GUI.
Where can I found it in order to check that the value is correct?

If you found something interesting in the listing, please let me know!
Thanks a lot.

JB
listing.txt (478 KB)

Hello,

The number of iterations might be limited by the fact the solver detects an divergence immediately, so stops there (you have a “nan” right away for the velocity).

Did you check the same computation without fan modeling ? I assume the error is due to the fan model, but it is safe to check anyways…

A debug build would be useful here, but otherwise, adding src/base_cs_fan.c to your user sources, could you add checks than elements are indeed selected for fans ?
In the “trunk” version (4.3 which will be released this summer), we added a visualisation of selected fan zones, but version 4.2 does not have this feature…

Regards,

Yvan

Hello,

I’ve run the simulation without fan modeling.
There is no problem. (it’s a pure internal fan in a cube, so here, I have a cube with wall BC on every faces and without any speed or moment in…). K and epsilon seems to converge correctly.

How can I add checks for selected cells for fan?
I’m really poor in C coding…

Thanks a lot

JB

Hello,

Some time ago I was working with other CFD software in a simulation with fan, so it is very interesting for me too.

I was following the comments of Yvan, I added the “cs_fan.c” to my 4.2.1 version and also have to add “cs_math.c” and “cs_math.h” because the functions cs_math_3_cross_product() and cs_math_3_norm() where not defined. Also I have to comment a pair of lines in cs_fan.c where is used

cell_fan_id

, may be this is a problem of my version.

/* Store the cell_fan_id in the postprocessing field */
  // cs_field_t *fans_ids = cs_field_by_name("fans_ids");
  //for (cs_lnum_t cell_id = 0; cell_id < n_ext_cells; cell_id++)
  // fans_ids->val[cell_id] = (cs_real_t)cell_fan_id[cell_id];

After that I have run a simple example ( A box with a fan inside) with your parameters. First of all I check your pressure curve and it has negative value of pressure, maybe there is a problem with the reference pressure. In general the curves are like this


I have run the example (no turbulence model) with your curve and with the used in the template of cs_user_parameters-fan.c without problem. I don’t understand very well the definition of “axial_torque” but it is a relation between the axial and radial velocities.

If you want I can upload the case

Regards.

Hello,

This seems interesting.

Just a note regarding the documentation (or lack thereof). The fan model was implemented and documented in an older code, and the model was recently added to Code_Saturne, but not documented yet. This documentation should be added to the code in the next few months, so tit should be either in version 4.3 (this summer), or 5.0 (next year, with a feature freeze around november).

The curves in our documentation are similar to the ones you posted, but I don’t have them here (only have access to the documentation system from the office, and have not made a local copy of this document), so I can’t tell right now if the representation is exactly the same, but it should ba at worse quite similar.

From memory, axial torque is an optional model parameter accounting for the fact that the rotation of the fans adds a rotation term to the flow (i.e. in a virtual “tube” containing the fan, the flow has a spiral, not straight movement).

Regards,

Yvan

Hello,

I suspected that the axial torque has the influence that you mention.

Wth respect to the pressure curve, is a function of volumetric flow rate or velocity, because I dont understand what means with flow in th doxygen doc.

Regards

Luciano

Thanks a lot Luciano for your interest in this post.

I’ll be very interested if you can post your case.
I’ve try to do the same copy as you, but it doesn’t works on my side (copying cs_mat.c, cs_math.h and cs_fan.c from the trunk repository to the SRC directory of the case, commenting same lines as you).

I’ve got a lot of compilation errors about the files above…

Where did you copy them?

Regards

JB

Ok,

The Monday morning dummy error…
I’ve copied the html page of the trunk instead of the cs_math.c and cs_math.h files…
Changed that and the case now compiles.

However, I still have a divergence in the first iteration, even without turbulence model.

Is your case a closed box? Or do you have Inlet and outlet on some walls?

Regards

JB

Hello,

In my test case I have set inlet and outlet.

I attach the case file, the mesh and the fan function. If this case run without problem post your case in order to help you to debug the problem.

Regards,

Luciano
Fan_test_case.tar.gz (1.23 MB)

Luciano,

Thanks a lot for your case

I’m gonna try it and let you know how it works on my CS.

Regards

JB

Hello Yvan,

Do you have any references for the implementation of the axial fan model? Because there are different alternatives depending on the blade, fan and hub radios. Also in the 2D modeling, there is weighting with the thickness and radios.

Thanks in advance.

Regards,
Luciano.

Hello Luciano,

I might be able to find an old report (in French) in which 1 or 2 pages describe the model and might provide a reference. I need to do a bit of archeology on my machine, but I’ll let you know.

Best regards,

Yvan

Hello Luciano,

Here is the theoretical documentation for Fan modeling in N3S, a predecessor to code_saturne, with the same fan models.
This is in French, but with diagrams and equations.

Best regards,

Yvan
fan-model-fr.pdf (883 KB)