Initial condition from fld or chkp

[vbaconnet]

Adds some functionalities outlined in #1283 for using fields as an initial condition.

Basically allows to initialize the simulation from fld or chkp files with the keyword case.<fluid/scalar>.initial_condition.type = field.

• For chkp the already existing interpolations (mesh to mesh and from one polynomial order to another) are usable.
• For fld we basically do the exact same interpolation workflow as chkp, but we wrap the findpts part in fld_data_file instead of having it in the reader. Chkp should perhaps be modified to do that as well if we decide to go that way.
• I am opening this as draft for others to add the missing functionalities if they so want, see the list below for the missing functionalities :)

Todo

So far this needs to be tested with:

• AMD GPUs
• Nvidia GPUs
• CPUs

Functionalities

• Read from chkp with or without interpolation
• Read from fld without interpolation
• Read from fld with interpolation (different polynomial order)
• Read from fld with interpolation (mesh to mesh)
• Read from fld with 2D extrusion (I think this is already quite a lot so probably we can skip this for now?)

Usage

chkp

We read the entire chkp file but only apply u,v,w,p and/or s.

"initial_condition": {
    "type": "field",
    "file_name": "fluid00001.chkp",
}

Any interpolation that chkp does when reading the file can also be done here, so essentially one can use a chkp file with mesh to mesh interpolation or just interpolation from different polynomial order (the latter doesn't need any additional keywords):

"initial_condition": {
    "type": "field",
    "file_name": "fluid00001.chkp",
    "previous_mesh": "mesh.nmsh",
    "tolerance": 1e-6
}

fld

Without interpolation

Here I see a two different use cases:

The first one would be where one has a series of field0.fld files and we select the index of the file we want to sample.

"initial_condition": {
    "type": "field",
    "file_name": "field0.fld" or "field0.nek5000",
    "sample_index": 5 <--- this will read field0.f00005
}

Without specifying sample_index we read the last fld file in the series by default, provided that a field0.nek5000 file exists, if not it will default to reading field0.f00000.

The second use case is if one only has a single field file for example myfield0.f00100, then one can directly use that file name and the sample_index will be retrieved from the extension f00100:

"initial_condition": {
    "type": "field",
    "file_name": "myfield0.f00100"
}

With interpolation

The interpolation between two different polynomial orders is done automatically without having to specify anything (same as for chkp).

In order to use the mesh to mesh interpolation, we use the interpolate keyword:

"initial_condition": {
    "type": "field",
    "file_name": "field0.fld" or "field0.nek5000",
    "sample_index": 5, <--- this will read field0.f00005
    "interpolate": true
}

When a fld series is given as file_name, we will by default read the x,y,z coordinates in the first fld file of the series. If we give a single .f***** file, it will try to look for the first file in the series (if it can find a .nek5000 file, otherwise it will read the f00000 file) for coordinates. If the coordinates are in a specific file we can also provide the index of the file with sample_mesh_index:

"initial_condition": {
    "type": "field",
    "file_name": "myfield0.f00100",
    "interpolate": true,
    "sample_mesh_index": 99 <--- means that the coordinates are in myfield0.f00099
}

Other changes

I added some logging for the initial condition so it becomes cleaner when using the fld reader that outputs some messages, that gives something like:

   ---Fluid initial condition----  
   Type: field
   File name: results/field0.fld
   Reading meta file for fld series
   Name: field0
   Start counter    :       0
   Number of samples:      17
   Reading fld file results/field0.f00000
   Reading fld file results/field0.f00016
  
   ---Scalar initial condition---  
   Type: field
   File name: results/field0.fld
   Sample index: 00005
   Reading meta file for fld series
   Name: field0
   Start counter    :       0
   Number of samples:      17
   Reading fld file results/field0.f00005

But it will also output all the values provided if one uses the other initial condition types like blasius, uniform, point_zone etc.

   ---Fluid initial condition----  
   Type: blasius
   delta       :   0.001000
   Approximation: sin
   Value: [  1.000000,  0.000000,  0.000000]

[vbaconnet]

Could we also init with chkp_data%init(s,s,s,s)? Or is it dangerous?

[Shiyu-Sandy-Du]

I don't think it would be dangerous anyway since in chkp_t%init(), all arguments are declared to be intent(in)

[vbaconnet]

Yes, that is true. But then chkp_t has pointers to u,v,w,p in checkpoint::chkp_init and the reader seems to write over them in chkp_file::chkp_file_read:

    call this%read_field(fh, byte_offset, u%x, nel)
! ...
    call this%read_field(fh, byte_offset, v%x, nel)
! ...
    call this%read_field(fh, byte_offset, w%x, nel)
! ...
    call this%read_field(fh, byte_offset, p%x, nel)
! ...
    if (read_scalar) then
       call this%read_field(fh, byte_offset, s%x, nel)
    end if

So I wonder if this means that a call to chkp_file%read(chkp_data) will overwrite u,v,w,p and s, potentially messing things up.

[adperezm]

I think this is a good place where maybe @MartinKarp knows from the top of his head.

[vbaconnet]

Actually now that I am testing just initializing the scalar from a checkpoint I see that it's not reading the scalar properly anymore, but u instead. Even with init(s,s,s,s), kind of weird

[timofeymukha]

Yeah, but the whole problem of dealing with raw binary data, where you have to care about the order, carefully compute offsets etc, will disappear. Adding something new to the checkpoint will also no longer be a big deal. Simcomps will just be able to drop in any data they may need for restart as hdf5 datasets.

[vbaconnet]

Is there something fundamentally different in what you are doing here?

Yes, I do not initialize any lag arrays in the chkp_t object before reading the chkp file, because I do not have access to those from flow_ic or scalar_ic, so chkp_file_read detects that there are lag arrays written in the file but since they are not associated from the chkp file it will skip the reading. When we restart from a chkp file it will read all of those so it does not skip anything

[vbaconnet]

For a normal restart you add the lag arrays in fluid_pnpn:

call this%chkp%add_lag(this%ulag, this%vlag, this%wlag)

but those are in fluid, which is not visible from scalar_ic

[MartinKarp]

Is the ic from chkp working now? Otherwise, as I wrote in my review I would perhaps skip it. I am also not entirely sure when one wants to set the ic from a checkpoint.

[vbaconnet]

Yes it's working but the size of lag arrays is hardcoded. Kind of annoying but that's the only way I can see without just skipping chkp entirely. Perhaps we don't even want that feature :)

[adperezm]

Great work!

A comment before reading the code:
I think it is very good also that one can specify directly the file name (avoiding .nek5000) files. In that sense, it could also be good that one gives that option for the file with the mesh, rather than just the index.

I think I see your reasoning, however, since you would restart in a folder that already has an up to date .nek5000 and all that. But it is good to have options.

[vbaconnet]

it could also be good that one gives that option for the file with the mesh, rather than just the index.

Actually yes it would make sense to be able to do it via a file too. Also I am now wondering if it even makes sense to have the option to give the index?

[adperezm]

I did not check the ic files since I am not familiar with them, but everything is very nice! just some minor comments.

[MartinKarp]

Good work! The only thing I would like to raise is whether we should support set ic from chkp. For the fld file I think it is great, but why would one not restart from the chkp?

Is there some interpolation we do when setting it as ic we don't do when doing the restart?

I would personally skip support for chkp, but now you have done a lot of work on it so.

Sorry for the late comment, but Im on vacation until the 11th.

[MartinKarp]

Yes, this is a problem I guess. The init from a dofmap is not needed. I think perhaps this would be better solved by adding another constructor for global_inteprolation that takes three arrays instead of a dofmap. This could be arranged similar to how the init for field is done if I remember correctly where we have two inits.

[vbaconnet]

why would one not restart from the chkp?

I guess it was just to have that option. In case the restart somehow does not work for any reason one can still use the fields for a first order restart. Thinking about some simulation component that does not restart properly for example. But in practice as you say this is probably quite rare.