16#include <Champ_Fonc_MED.h>
21#include <medcoupling++.h>
22#include <TRUST_2_MED.h>
23#include <Comm_Group_MPI.h>
25#include <MEDLoader.hxx>
26#include <MEDCouplingField.hxx>
27#include <MEDCouplingFieldDouble.hxx>
29#include <ParaMEDFileMesh.hxx>
31#pragma GCC diagnostic ignored "-Wreorder"
32#include <MEDFileField.hxx>
33#include <communications.h>
35using MEDCoupling::MEDCouplingField;
36using MEDCoupling::MEDCouplingFieldDouble;
37using MEDCoupling::MCAuto;
38using MEDCoupling::DataArrayIdType;
39using MEDCoupling::GetTimeAttachedOnFieldIteration;
40using MEDCoupling::GetAllFieldNamesOnMesh;
41using MEDCoupling::MEDFileFieldMultiTS;
42using MEDCoupling::MEDFileField1TS;
43using MEDCoupling::MEDFileMesh;
86 bool nom_decoup_lu =
false;
90 if (chaine_lue ==
"{")
98 if (motlu ==
"{") cnt++;
99 if (motlu ==
"}") cnt --;
100 s +=
Nom(
" ") + motlu;
107 constexpr double MIN_INF = -1.0e+300;
109 param.lire_avec_accolades(is2);
116 Cerr <<
"ERROR: reading 'Champ_fonc_MED': no time was provided! Either 'last_time' or 'time' parameter must be specified!" << finl;
121 Cerr <<
"ERROR: reading 'Champ_fonc_MED': both 'time' and 'last_time' were provided! Choose only one!" << finl;
128 Cerr <<
"ERROR: reading 'Champ_fonc_MED': Expected opening brace '{' - are you using the new syntax?" << finl;
129 Cerr <<
"If you are still using the old syntax (before TRUST v1.9.3), \nyou can use -convert_data option of your application script:" << finl;
143 int multiple_med = 0;
145 tmp.prefix(
"0001.med");
158 multiple_med =
mp_max(multiple_med);
177 Cerr <<
"Warning, you can't use use_existing_domain on a partitionned domain like " <<
nom_dom_ << finl;
178 Cerr <<
"It is not parallelized yet... So we use MED mesh, which is not optimal." << finl;
179 Cerr <<
"Try suppress use_existing_domain option." << finl;
187 Cerr<<
nom_dom_<<
" was not read into the med file because it already exists." <<finl;
196 Cerr <<
"Checking whether domain in the file "<<
nom_fichier_med_<<
" and domain "<<
nom_dom_<<
" are the same (coords,connectivity) ..."<<finl;
200 liremed.retrieve_MC_objects();
201 const MEDCouplingUMesh* new_um = liremed.get_mc_mesh();
202 const MEDCouplingUMesh* root_um = le_domaine.get_mc_mesh();
205 DataArrayIdType *dnup1=
nullptr, *dnup2=
nullptr;
213 if (dnup2) dnup2->decrRef();
214 if (dnup1) dnup1->decrRef();
216 catch(INTERP_KERNEL::Exception& e)
218 Cerr <<
"Comparison of the two domains failed, with the following message from MEDCoupling::checkGeoEquivalWith()" << finl;
219 Cerr << e.what() << finl;
230 Cerr<<
"INFO: You can toggle the flag 'use_existing_domain' in 'Champ_Fonc_MED' to optimize reading since it seems that the domain already exists."<<finl;
236 liremed.lire_geom(
false);
240 DoubleTab BB1 =
dom_med_.getBoundingBox();
244 for (
int j=0; j<2; j++)
245 if (est_different(BB1(i,j), BB2(i,j)))
247 Cerr <<
"Error, Champ_Fonc_MED can't read the partitionned MED files." << finl;
248 Cerr <<
"=> it seems that the MED partition is different than the '" << dom_calcul.
le_nom() <<
"' domain partition." << finl;
251 Cerr <<
"Ok MED partition matches the domain partition so reading multiple MED files..." << finl;
260 Cout <<
"The resumption time is "<<
temps_<<finl;
277 for (
int i = 0; i < dec_size; i++)
280 else if(
loc_ ==
"som")
284 std::vector<std::set<int>> dec(nbNodes);
285 for(
int n=0; n<nbNodes; n++)
288 fdec >> node >> size;
289 for(
int p=0; p<size; p++)
293 dec[node].insert(proc);
296 dec_size = (int)dec.size();
297 for (
int n=0; n < dec_size; n++)
298 for(std::set<int>::iterator it = dec[n].begin(); it!=dec[n].
end(); ++it)
303 Cerr <<
"Champ_Fonc_MED on parallel domain : decoup file only handled for field located on nodes or on cells!" << finl;
311 for (
int i=0; i<nb_item; i++)
312 filter.push_back(first_item);
314 Cerr <<
"Creating a filter to access efficiently values in " <<
nom_fichier_med_ << finl;
317 Cerr <<
"Champ_Fonc_MED on parallel domain : inconsistency between filter and domain (not the same number of entities)!" << finl;
323 Cerr <<
"Champ_Fonc_MED on existing domain : inconsistency between domain file and field (not the same number of entities)!" << finl;
339 int given_it = -(int) t - 1;
347med_geometry_type type_geo_trio_to_type_med(
const Nom& type);
349INTERP_KERNEL::NormalizedCellType type_geo_trio_to_type_medcoupling(
const Nom& type,
int& mesh_dimension);
359 le_champ().Champ_Fonc_base::mettre_a_jour(t);
369 bool search_field =
true;
383 Cout <<
"We assume that the field " << fieldName <<
" is stationary." << finl;
384 search_field =
false;
390 temps_sauv_ = lire_temps_champ(fileName, fieldName);
394 for (given_it = 0; given_it < (int) nn; given_it++)
396 if (given_it == (
int)nn)
398 Cerr <<
"Error. Time " << t <<
" not found in the times list of the " << fileName <<
" file" << finl;
399 for (
unsigned int n=0; n<nn; n++)
404 int iteration = time_steps_[given_it].first;
405 int order = time_steps_[given_it].second;
410 MEDFileUtilities::AutoFid fid(MEDCoupling::OpenMEDFileForRead(fileName));
411 int mesh_dimension = -1;
412 std::vector<std::pair<MEDCoupling::TypeOfField,INTERP_KERNEL::NormalizedCellType>> tmp= { { field_type, type_geo_trio_to_type_medcoupling(mon_dom->type_elem()->que_suis_je(), mesh_dimension) } };
413 INTERP_KERNEL::AutoCppPtr<MEDCoupling::MEDFileEntities> entities(MEDCoupling::MEDFileEntities::BuildFrom(&tmp));
414 MCAuto<MEDCoupling::MEDFileAnyTypeField1TS> fieldFile = MEDCoupling::MEDFileAnyTypeField1TS::NewAdv(fid, fieldName, iteration, order, entities,
filter);
415 MCAuto<MEDCoupling::MEDFileField1TS> ret(MEDCoupling::DynamicCast<MEDCoupling::MEDFileAnyTypeField1TS, MEDCoupling::MEDFileField1TS>(fieldFile));
418 MEDCoupling::DataArrayDouble *field_values = ret->getUndergroundDataArray();
419 std::copy(field_values->begin(),field_values->end(),
425 MCAuto<MEDCouplingField> ffield;
426 bool has_field =
ffield_!=
nullptr;
427 if (!has_field) ffield = lire_champ(fileName, meshName, fieldName, iteration, order);
428 Cerr <<
" at time " << t <<
" ... " << finl;
429 MEDCouplingFieldDouble * field =
dynamic_cast<MEDCouplingFieldDouble *
>((MEDCouplingField *)(has_field?
ffield_:ffield));
432 Cerr <<
"ERROR reading MED field! Not a MEDCouplingFieldDouble!!" << finl;
436 assert((
int) field->getNumberOfComponents() ==
438 MEDCoupling::DataArrayDouble *field_values = field->getArray();
439 std::copy(field_values->begin(),field_values->end(),
446 le_champ().Champ_Fonc_base::mettre_a_jour(t);
462 int mesh_dimension = -1;
463 cell_type = type_geo_trio_to_type_medcoupling(type_elem, mesh_dimension);
464 if (localisation !=
Nom())
466 if (localisation ==
"som")
468 field_type = MEDCoupling::ON_NODES;
469 if ((cell_type == INTERP_KERNEL::NORM_QUAD4) || (cell_type == INTERP_KERNEL::NORM_HEXA8))
470 type_champ =
"Champ_Fonc_Q1_MED";
472 type_champ =
"Champ_Fonc_P1_MED";
474 else if (localisation ==
"elem")
476 field_type = MEDCoupling::ON_CELLS;
477 type_champ =
"Champ_Fonc_P0_MED";
481 Cerr << localisation <<
" localization unknown." << finl;
486 std::string fileName = nom_fic.
getString();
488 Noms fieldNamesGuess;
490 if (localisation !=
Nom())
496 std::vector<std::string> fieldNames = GetAllFieldNamesOnMesh(fileName, meshName);
497 for (
int i=0; i<fieldNamesGuess.size(); i++)
499 if (std::find(fieldNames.begin(), fieldNames.end(), fieldNamesGuess[i].getString()) != fieldNames.end())
508 Cerr <<
"Unable to find into file '" << fileName <<
"' a field named like :" << finl;
509 Cerr <<
" " << fieldNamesGuess << finl;
510 Cerr <<
"supported by mesh '" << meshName <<
"'" << finl;
512 Cerr <<
"This file contains the field(s) named:" << finl;
513 for (
unsigned i=0; i<fieldNames.size(); i++)
514 Cerr << fieldNames[i] << finl;
522 lire_donnees_champ(fileName,meshName,fieldName,temps_sauv,size,nbcomp,type_champ);
525 vrai_champ_.typer(type_champ);
543ArrOfDouble Champ_Fonc_MED::lire_temps_champ(
const std::string& fileName,
const std::string& fieldName)
545 ArrOfDouble temps_sauv;
546 MCAuto<MEDFileFieldMultiTS> ft1(MEDFileFieldMultiTS::New(fileName, fieldName));
547 std::vector<double> tps;
548 time_steps_ = ft1->getTimeSteps(tps);
549 unsigned int nn = (unsigned)tps.size();
551 for (
unsigned it = 0; it < nn; it++)
552 temps_sauv[it] = tps[it];
561void Champ_Fonc_MED::lire_donnees_champ(
const std::string& fileName,
const std::string& meshName,
const std::string& fieldName,
562 ArrOfDouble& temps_sauv,
int& size,
int& nbcomp,
Nom& type_champ)
564 temps_sauv = lire_temps_champ(fileName, fieldName);
566 int last_iter = time_steps_[nn-1].first;
567 int last_order = time_steps_[nn-1].second;
571 MCAuto<MEDFileField1TS> field = MEDFileField1TS::New(fileName, fieldName, last_iter, last_order);
572 nbcomp = (int) field->getNumberOfComponents();
577 MCAuto<MEDFileMesh> mesh = MEDFileMesh::New(fileName, meshName, field->getMeshIteration(), field->getMeshOrder());
578 size = field_type == MEDCoupling::ON_NODES ? (int)mesh->getNumberOfNodes() : (int)mesh->getNumberOfCellsAtLevel(0);
583 ffield_ = lire_champ(fileName, meshName, fieldName, last_iter, last_order);
584 MEDCouplingFieldDouble * field =
dynamic_cast<MEDCouplingFieldDouble *
>((MEDCouplingField *)
ffield_);
587 Cerr <<
"ERROR reading MED field! Not a MEDCouplingFieldDouble!!" << finl;
590 size = (int)field->getNumberOfTuplesExpected();
591 nbcomp = (int)field->getNumberOfComponents();
595 if (field_type == MEDCoupling::ON_CELLS)
596 type_champ =
"Champ_Fonc_P0_MED";
597 else if (field_type == MEDCoupling::ON_NODES)
598 type_champ = cell_type == INTERP_KERNEL::NORM_QUAD4 || cell_type == INTERP_KERNEL::NORM_HEXA8 ?
"Champ_Fonc_Q1_MED" :
"Champ_Fonc_P1_MED";
618MCAuto<MEDCouplingField> Champ_Fonc_MED::lire_champ(
const std::string& fileName,
const std::string& meshName,
619 const std::string& fieldName,
const int iteration,
const int order)
624 Cerr <<
"Reading" << (fast ?
" (fast)" :
"") <<
" the field " << fieldName <<
" on the " << meshName <<
" mesh into " << fileName <<
" file" << finl;
625 MCAuto<MEDCouplingField> ffield;
628 MCAuto<MEDFileField1TS> file = MEDFileField1TS::New(fileName, fieldName, iteration, order);
629 ffield = file->getFieldOnMeshAtLevel(field_type,
domaine().get_mc_mesh(), 0);
633 ffield = ReadField(field_type, fileName, meshName, 0, fieldName, iteration, order);
DoubleTab & valeurs() override
Surcharge Champ_base::valeurs() Renvoie le tableau des valeurs.
classe Champ_Fonc_MED Load a field from a MED file for a given time.
const Domaine & domaine() const override
std::vector< trustIdType > filter
virtual void lire(double tps, int given_iteration=-1)
MCAuto< MEDCouplingField > ffield_
DoubleTab & valeurs() override
Surcharge Champ_base::valeurs() Renvoie le tableau des valeurs.
int creer(const Nom &, const Domaine &dom, const Motcle &localisation, ArrOfDouble &temps_sauv)
void mettre_a_jour(double) override
Mise a jour en temps du champ.
Domaine_VF_inst domainebidon_inst
virtual const Champ_Fonc_base & le_champ() const
virtual void set_param(Param ¶m) const override
Nom nom_champ_dans_fichier_med_
bool use_existing_domain_
classe Champ_Fonc_base Classe de base des champs qui sont fonction d'une grandeur calculee
void associer_domaine_dis_base(const Domaine_dis_base &) override
void mettre_a_jour(double temps) override
Mise a jour en temps du champ.
int fixer_nb_valeurs_nodales(int nb_noeuds) override
Fixe le nombre de degres de liberte par composante.
virtual double changer_temps(const double t)
Fixe le temps auquel se situe le champ.
void corriger_unite_nom_compo()
cette methode va fixer les unites et le nom des compos elle n'est pas const en realite !...
DoubleTab getBoundingBox() const
void build_mc_mesh(bool virt=false) const
Build the MEDCoupling mesh corresponding to the TRUST mesh.
bool is_mc_mesh_ready() const
int_t nb_som() const
Renvoie le nombre de sommets du domaine.
const Nom & le_nom() const override
Donne le nom de l'Objet_U Methode a surcharger : renvoie "neant" dans cette implementation.
virtual int ouvrir(const char *name, IOS_OPEN_MODE mode=ios::in)
Class defining operators and methods for all reading operation in an input flow (file,...
virtual void fixer_nb_comp(int i)
Fixe le nombre de composantes du champ.
const Nom & le_nom() const override
Renvoie le nom du champ.
void nommer(const Nom &) override
Donne un nom au champ.
static int objet_existant(const Nom &)
Renvoie 1 si l'objet existe, 0 sinon voir Interprete_bloc::objet_global_existant().
Une chaine de caractere (Nom) en majuscules.
class Nom Une chaine de caractere pour nommer les objets de TRUST
const std::string & getString() const
const Nom & le_nom() const override
Renvoie *this;.
Un tableau de chaine de caracteres (VECT(Nom)).
int contient_(const char *const ch) const
const Interprete & interprete() const
const Nom & que_suis_je() const
renvoie la chaine identifiant la classe.
virtual Entree & readOn(Entree &)
Lecture d'un Objet_U sur un flot d'entree Methode a surcharger.
static double precision_geom
static const Nom & nom_du_cas()
Renvoie une reference constante vers le nom du cas.
virtual Sortie & printOn(Sortie &) const
Ecriture de l'objet sur un flot de sortie Methode a surcharger.
Helper class to factorize the readOn method of Objet_U classes.
void ajouter_flag(const char *keyword, const bool *value)
Register a boolean flag whose mere presence switches it to true.
void ajouter(const char *keyword, const int *value, Param::Nature nat=Param::OPTIONAL)
Register an integer parameter.
static trustIdType mppartial_sum(trustIdType i)
Calul de la somme partielle de i sur les processeurs 0 a me()-1 (renvoie 0 sur le processeur 0).
static double mp_max(double)
static bool is_parallel()
static int nproc()
renvoie le nombre de processeurs dans le groupe courant Voir Comm_Group::nproc() et PE_Groups::curren...
static int me()
renvoie mon rang dans le groupe de communication courant.
static void exit(int exit_code=-1)
Routine de sortie de TRUST dans une region Kokkos.
Classe de base des flux de sortie.
_SIZE_ size_array() const
void resize_array(_SIZE_ new_size, RESIZE_OPTIONS opt=RESIZE_OPTIONS::COPY_INIT)
_SIZE_ dimension(int d) const