1 # -*- coding: utf-8 -*-
3 # Copyright (C) 2016-2024 EDF R&D
5 # This library is free software; you can redistribute it and/or
6 # modify it under the terms of the GNU Lesser General Public
7 # License as published by the Free Software Foundation; either
8 # version 2.1 of the License.
10 # This library is distributed in the hope that it will be useful,
11 # but WITHOUT ANY WARRANTY; without even the implied warranty of
12 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13 # Lesser General Public License for more details.
15 # You should have received a copy of the GNU Lesser General Public
16 # License along with this library; if not, write to the Free Software
17 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 # Author: Jean-Philippe Argaud, jean-philippe.argaud@edf.fr, EDF R&D
20 # Author: Luis Corona Mesa-Molles, luis.corona-mesa-molles@edf.fr, EDF R&D
22 __all__ = ["Calibration", "name", "version", "year"]
24 # ==============================================================================
26 name = "Modelica and Dymola Calibration Tools"
27 version = "1.0.9.7.0" # "x.x"+"adao version"
30 # ==============================================================================
31 # Default configuration
32 # ---------------------
33 import os, sys, io, shutil, numpy, scipy, tempfile, time, logging, pandas, subprocess, copy, csv
34 from datetime import timedelta
35 from datetime import datetime
39 from adao import adaoBuilder
40 from daCore.Interfaces import ImportFromFile, ImportScalarLinesFromFile
42 raise ImportError("ADAO module not found, please install it first.")
45 from buildingspy.io.outputfile import Reader
47 raise ImportError("buildingspy module not found, please install it first.")
50 from pst4mod.modelica_libraries.automatic_simulation import Automatic_Simulation
51 from pst4mod.modelica_libraries.around_simulation import Around_Simulation
52 from pst4mod.modelica_libraries import write_in_dsin
53 from pst4mod.modelica_libraries.change_working_point import Dict_Var_To_Fix
54 from pst4mod.modelica_libraries.change_working_point import Working_Point_Modification
56 raise ImportError("modelica_libraries module not found, please install it first.")
59 from fmpy import simulate_fmu
60 from fmpy.fmi1 import FMICallException
62 # raise ImportError("fmpy library not found, please install it first")
67 "DirectoryModels" : "Models",
68 "DirectoryMeasures" : "Measures",
69 "DirectoryMethods" : "Methods",
70 "DirectoryResults" : "Results",
71 "Launcher" : "configuration.py",
74 __paths.append( os.path.abspath(os.path.join(os.path.dirname(__file__), '..','..','buildingspy-1.5.0')) )
75 __paths.append( os.path.abspath(os.path.join(os.path.dirname(__file__), '..','..','Packages')) )
76 for __path in __paths:
77 if os.path.exists(__path): sys.path.append(__path)
79 # ==============================================================================
80 class Calibration(object):
81 "ADAO case based parameters inverse estimation tools"
82 def __init__(self, Name = "Calibration case", SaveStdoutIn = None, Verbose = False):
84 Name : name as a string
85 SaveAllStdoutIn: filename to be used for stdout stream
86 Verbose : verbose output
88 self.__name = str(Name)
94 self._setConfigurationDefaults()
95 self.__verbose = bool(Verbose)
96 self.__stdoutid = sys.stdout
97 if SaveStdoutIn is not None:
98 sys.stdout = open(SaveStdoutIn,"w")
100 __msg = "[VERBOSE] %s"%self.__name
103 print(" %s"%("-"*len(__msg),))
105 VersionsLogicielles()
107 def _setConfigurationDefaults(self, Configuration = None):
109 Impose case directory and files structure configuration defaults based
110 on argument dictionary
112 if Configuration is None or type(Configuration) is not dict:
113 Configuration = _configuration
114 for __k,__v in Configuration.items():
115 setattr(self, __k, __v)
117 def _setModelTmpDir(self, __model_name="dsin.txt", __model_format="DYMOSIM", __model_dest = "dsin.txt", __before_tmp = None, __suffix=None):
118 if __model_name is None:
119 raise ValueError('Model file or directory has to be set and can not be None')
120 elif os.path.isfile(__model_name):
121 __model_nam = os.path.basename(__model_name)
122 __model_dir = os.path.abspath(os.path.dirname(__model_name))
123 __model_dst = __model_dest
124 elif os.path.isdir(__model_name):
125 __model_nam = "dsin.txt"
126 __model_dir = os.path.abspath(__model_name)
127 __model_dst = __model_dest
129 raise ValueError('Model file or directory not found using %s'%str(__model_name))
131 if __before_tmp is not None: # Mettre "../.." si nécessaire
132 __mtmp = os.path.join(__model_dir, __before_tmp, "tmp")
134 __mtmp = os.path.join(__model_dir, "tmp")
135 if not os.path.exists(__mtmp):
137 __prefix = time.strftime('%Y%m%d_%Hh%Mm%Ss_tmp_',time.localtime())
138 if __suffix is not None:
139 __prefix = __prefix + __suffix + "_"
140 __ltmp = tempfile.mkdtemp( prefix=__prefix, dir=__mtmp )
142 for sim, pmf, dst in (
143 (__model_nam,__model_format,__model_dst),
144 ("dymosim.exe","DYMOSIM","dymosim.exe"),
145 ("pythonsim.exe","PYSIM","pythonsim.exe"),
147 # Recherche un fichier de données ou un simulateur dans cet ordre
149 if os.path.isfile(os.path.join(__model_dir, sim)) and (pmf.upper() in [__model_format, "GUESS"]):
151 os.path.join(__model_dir, sim),
152 os.path.join(__ltmp, dst)
155 # _secure_copy_textfile(
156 # os.path.join(__model_dir, sim),
157 # os.path.join(__ltmp, dst)
163 def _setModelTmpDir_simple(self, __model_name="dsin.txt", __model_format="DYMOSIM", __model_dest = "dsin.txt"):
164 if __model_name is None:
165 raise ValueError('Model file or directory has to be set and can not be None')
166 elif os.path.isfile(__model_name):
167 if __model_format == "OPENMODELICA": #same structure as in the directory case
168 __model_dir = os.path.abspath(os.path.dirname(__model_name))
170 __model_nam_1 = os.path.basename(__model_name)
171 __model_nam_2 = __model_nam_1[:-9]+'.exe'
172 __model_nam = [__model_nam_1, __model_nam_2] #the first found model is kept
174 __model_nam_3 = __model_nam_1[:-9]+'_info.json' #check if this file exists as well
175 if os.path.exists(os.path.join(__model_dir,__model_nam_3)):
176 __model_nam.append(__model_nam_3) #get the three files necessar for simulation
177 __model_dst = __model_nam #the same file name is kept
180 __model_nam = os.path.basename(__model_name)
181 __model_dir = os.path.abspath(os.path.dirname(__model_name))
182 __model_dst = os.path.basename(__model_dest)
183 elif os.path.isdir(__model_name):
184 if __model_format == "DYMOSIM":
185 __model_nam = "dsin.txt"
186 __model_dir = os.path.abspath(__model_name)
187 __model_dst = os.path.basename(__model_dest)
189 elif __model_format == "FMI" or __model_format == "FMU":
190 __model_dir = os.path.abspath(__model_name)
191 __model_nam = [files for files in os.listdir(__model_dir) if files[-4:] =='.fmu'][0] #the first found model is kept
192 __model_dst = __model_nam #the same file name is kept
194 elif __model_format == "OPENMODELICA":
195 __model_dir = os.path.abspath(__model_name)
196 __model_nam_1 = [files for files in os.listdir(__model_dir) if files[-4:] =='.xml'][0] #the first found model is kept, .xml extension is considered
197 __model_nam_2 = __model_nam_1[:-9]+'.exe'
198 __model_nam = [__model_nam_1, __model_nam_2] #the first found model is kept
200 __model_nam_3 = __model_nam_1[:-9]+'_info.json' #check if this file exists as well
201 if os.path.exists(os.path.join(__model_dir,__model_nam_3)):
202 __model_nam.append(__model_nam_3) #get the three files necessar for simulation
203 __model_dst = __model_nam #the same file name is kept
206 raise ValueError('If the model is not provided in DYMOSIM, OpenModelica or FMI format, the name of the model must refer to a file name (in a relative way with respect to the configuration file)')
208 raise ValueError('Model file or directory not found using %s'%str(__model_name))
210 __mtmp = os.path.join(__model_dir, "tmp")
211 if not os.path.exists(__mtmp):
215 if type(__model_nam) == list: #it means that it is an OpenModelica model
216 for sim_element in __model_nam:
218 os.path.join(__model_dir, sim_element),
219 os.path.join(__ltmp, sim_element) #the same name is kept for the files
222 for sim, pmf, dst in (
223 (__model_nam,__model_format,__model_dst),
224 ("dsin.txt","DYMOSIM","dsin.txt"),
225 ("pythonsim.exe","PYSIM","pythonsim.exe"),
227 # Recherche un fichier de données ou un simulateur dans cet ordre
228 if os.path.isfile(os.path.join(__model_dir, sim)) and (pmf.upper() in [__model_format, "GUESS"]):
231 os.path.join(__model_dir, sim),
232 os.path.join(__ltmp, dst)
234 # _secure_copy_textfile(
235 # os.path.join(__model_dir, sim),
236 # os.path.join(__ltmp, dst)
242 def _get_Name_ModelTmpDir_simple(self, __model_name="dsin.txt"):
243 if __model_name is None:
244 raise ValueError('Model file or directory has to be set and can not be None')
245 elif os.path.isfile(__model_name):
246 __model_dir = os.path.abspath(os.path.dirname(__model_name))
247 elif os.path.isdir(__model_name):
248 __model_dir = os.path.abspath(__model_name)
250 raise ValueError('Model file or directory not found using %s'%str(__model_name))
252 __mtmp = os.path.join(__model_dir, "tmp")
256 def _setModelTmpDir_REF_CALCULATION(self, __model_name="dsin.txt", __model_format="DYMOSIM", __suffix=None, __model_dest = "dsin.txt"):
257 if __model_name is None:
258 raise ValueError('Model file or directory has to be set and can not be None')
259 elif os.path.isfile(__model_name):
260 __model_dir = os.path.abspath(os.path.dirname(__model_name))
261 __model_nam = os.path.basename(__model_name)
262 elif os.path.isdir(__model_name):
263 __model_dir = os.path.abspath(__model_name)
264 __model_nam = "dsin.txt"
266 raise ValueError('Model file or directory not found using %s'%str(__model_name))
268 raise ValueError('Observation files must be named')
270 __mtmp = os.path.join(__model_dir, "tmp")
271 if not os.path.exists(__mtmp):
273 for sim, pmf in ((__model_nam,__model_format), ("dymosim.exe","DYMOSIM"), ("pythonsim.exe","PYSIM")):
274 # Recherche un simulateur dans cet ordre
275 if os.path.isfile(os.path.join(__model_dir, sim)) and (pmf.upper() in [__model_format, "GUESS"]):
276 __pref = time.strftime('%Y%m%d_%Hh%Mm%Ss_REF_CALCULATION_',time.localtime())
277 if __suffix is not None:
278 __pref = __pref+__suffix+"_"
279 __ltmp = tempfile.mkdtemp( prefix=__pref, dir=__mtmp )
281 os.path.join(__model_dir, sim),
282 os.path.join(__ltmp, __model_dest)
284 # _secure_copy_textfile(
285 # os.path.join(__model_dir, sim),
286 # os.path.join(__ltmp, __model_dest)
294 raise ValueError("Model simulator not found using \"%s\" path and \"%s\" format"%(__model_name,__model_format))
297 def describeMeasures(self, SourceName, Variables, Format="Guess"):
298 "To store description of measures"
299 self.__measures["SourceName"] = str(SourceName) # File or data name
300 self.__measures["Variables"] = tuple(Variables) # Name of variables
301 self.__measures["Format"] = str(Format) # File or data format
302 if self.__verbose: print(" [VERBOSE] Measures described")
303 return self.__measures
305 def describeMultiMeasures(self, SourceNames, Variables, Format="Guess"):
306 "To store description of measures"
307 self.__measures["SourceName"] = tuple(SourceNames) # List of file or data name
308 self.__measures["Variables"] = tuple(Variables) # Name of variables
309 self.__measures["Format"] = str(Format) # File or data format
310 if self.__verbose: print(" [VERBOSE] Measures described")
311 return self.__measures
313 def describeModel(self, SourceName, VariablesToCalibrate, OutputVariables, Format="Guess"):
314 "To store description of model"
315 self.__model["SourceName"] = str(SourceName) # Model command file
316 self.__model["VariablesToCalibrate"] = tuple(VariablesToCalibrate) # Name of variables
317 self.__model["OutputVariables"] = tuple(OutputVariables) # Name of variables
318 self.__model["Format"] = str(Format) # File or model format
319 if os.path.isfile(SourceName):
320 self.__model["Verbose_printing_name"] = os.path.basename(SourceName)
321 elif os.path.isdir(SourceName):
322 if self.__model["Format"] == "DYMOSIM":
323 self.__model["Verbose_printing_name"] = "dsin.txt + dymosim.exe in " + os.path.basename(SourceName)
324 elif self.__model["Format"] == "FMI" or self.__model["Format"] == "FMU":
325 self.__model["Verbose_printing_name"] = [files for files in os.listdir(os.path.abspath(SourceName)) if files[-4:] =='.fmu'][0]
326 elif self.__model["Format"] == "OPENMODELICA":
327 xml_file_name = [files for files in os.listdir(os.path.abspath(SourceName)) if files[-4:] =='.xml'][0]
328 exe_file_name = xml_file_name[:-9] + '.exe'
329 json_file_name = xml_file_name[:-9] + '_info.json'
330 if os.path.exists(os.path.join(os.path.abspath(SourceName),json_file_name)):
331 self.__model["Verbose_printing_name"] = xml_file_name + ' ' + exe_file_name + ' ' + json_file_name
333 self.__model["Verbose_printing_name"] = xml_file_name + ' ' + exe_file_name + ' '
335 print(" [VERBOSE] Model described: ", "/!\\ Unknown ModelFormat /!\\ ")
337 if self.__verbose: print(" [VERBOSE] Model described: ", self.__model["Verbose_printing_name"] )
340 def describeLink(self, MeasuresNames, ModelName, LinkName, LinkVariable, Format="Guess"):
341 "To store description of link between measures and model input data"
342 self.__link["SourceName"] = str(LinkName) # File or link name
343 self.__link["MeasuresNames"] = tuple(MeasuresNames) # List of file or name
344 self.__link["ModelName"] = tuple(ModelName) # File or model name
345 self.__link["Variable"] = str(LinkVariable) # Name of variable names column
346 self.__link["Format"] = str(Format) # File or model format
347 if self.__verbose: print(" [VERBOSE] Link between model and multiple measures described")
350 def describeMethod(self, SourceName, MethodFormat="PY", BackgroundName=None, BackgroundFormat="Guess"):
351 "To store description of calibration method"
352 self.__method["SourceName"] = str(SourceName) # File or method name
353 self.__method["Format"] = str(MethodFormat) # File or model format
354 self.__method["BackgroundName"] = str(BackgroundName) # File or background name
355 self.__method["BackgroundFormat"] = str(BackgroundFormat) # File or model format
356 if self.__verbose: print(" [VERBOSE] Optimization method settings described")
359 def describeResult(self, ResultName, Level=None, Format="TXT"):
360 "To store description of results"
361 self.__result["SourceName"] = str(ResultName) # File or method name
362 self.__result["Level"] = str(Level) # Level = "Final", "IntermediaryFinal"
363 self.__result["Format"] = str(Format) # File or model format
364 if self.__verbose: print(" [VERBOSE] Results settings described")
368 ModelName = None, # Model command file
369 ModelFormat = "Guess",
370 DataName = None, # Measures
371 DataFormat = "Guess",
372 BackgroundName = None, # Background
373 BackgroundFormat = "Guess",
374 MethodName = None, # Assimilation
375 MethodFormat = "Guess",
376 VariablesToCalibrate = None,
377 OutputVariables = None,
378 ResultName = "results.txt", # Results
379 ResultFormat = "Guess",
384 General method to set and realize a calibration (optimal estimation
385 task) of model parameters using measures data.
388 self.__verbose = True
389 print(" [VERBOSE] Verbose mode activated")
391 self.__verbose = False
392 if DataName is not None:
393 # On force l'utilisateur à nommer dans son fichier de mesures
394 # les variables avec le même nom que des sorties dans le modèle
395 self.describeMeasures(DataName, OutputVariables, DataFormat)
396 if ModelName is not None:
397 self.describeModel(ModelName, VariablesToCalibrate, OutputVariables, ModelFormat)
398 if MethodName is not None:
399 self.describeMethod(MethodName, MethodFormat, BackgroundName, BackgroundFormat)
400 if ResultName is not None:
401 self.describeResult(ResultName, ResultLevel, ResultFormat)
403 ONames, Observations = _readMeasures(
404 self.__measures["SourceName"],
405 self.__measures["Variables"],
406 self.__measures["Format"],
408 Algo, Params, CovB, CovR = _readMethod(
409 self.__method["SourceName"],
410 self.__method["Format"],
413 if self.__method["BackgroundFormat"] == "DSIN": __bgfile = self.__model["SourceName"]
414 if self.__method["BackgroundFormat"] == "ADAO": __bgfile = self.__method["SourceName"]
415 if self.__method["BackgroundFormat"] == "USER": __bgfile = self.__method["BackgroundName"]
416 if self.__method["BackgroundFormat"] == "Guess":
417 if self.__method["BackgroundName"] is not None: __bgfile = self.__method["BackgroundName"]
418 if self.__method["SourceName"] is not None: __bgfile = self.__method["SourceName"]
419 if self.__model["SourceName"] is not None: __bgfile = self.__model["SourceName"]
420 BNames, Background, Bounds = _readBackground(
422 self.__model["VariablesToCalibrate"],
423 self.__method["BackgroundFormat"],
425 if "Bounds" not in Params: # On force la priorité aux bornes utilisateur
426 Params["Bounds"] = Bounds
428 print(" [VERBOSE] Measures read")
429 print(" [VERBOSE] Optimization information read")
430 if "MaximumNumberOfSteps" in Params:
431 print(" [VERBOSE] Maximum possible number of iteration:",Params['MaximumNumberOfSteps'])
432 print(" [VERBOSE] Background read:",Background)
433 print(" [VERBOSE] Bounds read:",Params['Bounds'])
435 if "DifferentialIncrement" not in Params:
436 Params["DifferentialIncrement"] = 0.001
437 if "StoreSupplementaryCalculations" not in Params:
438 Params["StoreSupplementaryCalculations"] = ["SimulatedObservationAtOptimum",]
439 if "StoreSupplementaryCalculations" in Params and \
440 "SimulatedObservationAtOptimum" not in Params["StoreSupplementaryCalculations"]:
441 Params["StoreSupplementaryCalculations"].append("SimulatedObservationAtOptimum")
442 if self.__verbose and "StoreSupplementaryCalculations" in Params and \
443 "CurrentOptimum" not in Params["StoreSupplementaryCalculations"]:
444 Params["StoreSupplementaryCalculations"].append("CurrentOptimum")
445 if self.__verbose and "StoreSupplementaryCalculations" in Params and \
446 "CostFunctionJAtCurrentOptimum" not in Params["StoreSupplementaryCalculations"]:
447 Params["StoreSupplementaryCalculations"].append("CostFunctionJAtCurrentOptimum")
449 def exedymosim( x_values_matrix ):
450 "Appel du modèle et restitution des résultats"
451 from buildingspy.io.outputfile import Reader
452 from pst4mod.modelica_libraries.automatic_simulation import Automatic_Simulation
454 # x_values = list(numpy.ravel(x_values_matrix) * Background)
455 x_values = list(numpy.ravel(x_values_matrix))
456 dict_inputs=dict(zip(VariablesToCalibrate,x_values))
457 # if Verbose: print(" Simulation for %s"%numpy.ravel(x_values))
459 simudir = self._setModelTmpDir(ModelName, ModelFormat)
460 auto_simul = Automatic_Simulation(
462 dict_inputs=dict_inputs,
465 without_modelicares=True) #linux=true is removed
466 auto_simul.single_simulation()
467 reader = Reader(os.path.join(simudir,'dsres.mat'),'dymola')
468 y_values = [reader.values(y_name)[1][-1] for y_name in OutputVariables]
469 y_values_matrix = numpy.asarray(y_values)
471 shutil.rmtree(simudir, ignore_errors=True)
472 return y_values_matrix
474 def exepython( x_values_matrix ):
475 "Appel du modèle et restitution des résultats"
476 # x_values = list(numpy.ravel(x_values_matrix))
477 # dict_inputs=dict(zip(VariablesToCalibrate,x_values))
479 simudir = self._setModelTmpDir(ModelName, ModelFormat)
480 auto_simul = Python_Simulation(
482 val_inputs=x_values_matrix,
485 verbose=self.__verbose)
486 y_values_matrix = auto_simul.single_simulation()
488 shutil.rmtree(simudir, ignore_errors=True)
489 return y_values_matrix
491 if self.__model["Format"].upper() in ["DYMOSIM", "GUESS"]:
492 simulateur = exedymosim
493 elif self.__model["Format"].upper() == "PYSIM":
494 simulateur = exepython
496 raise ValueError("Unknown model format \"%s\""%self.__model["Format"].upper())
498 __adaocase = adaoBuilder.New()
499 __adaocase.set( 'AlgorithmParameters', Algorithm=Algo, Parameters=Params )
500 __adaocase.set( 'Background', Vector=Background)
501 __adaocase.set( 'Observation', Vector=Observations )
502 if type(CovB) is float:
503 __adaocase.set( 'BackgroundError', ScalarSparseMatrix=CovB )
505 __adaocase.set( 'BackgroundError', Matrix=CovB )
506 if type(CovR) is float:
507 __adaocase.set( 'ObservationError', ScalarSparseMatrix=CovR )
509 __adaocase.set( 'ObservationError', Matrix=CovR )
510 __adaocase.set( 'ObservationOperator', OneFunction=simulateur, Parameters = {"DifferentialIncrement":Params["DifferentialIncrement"]} )
512 __adaocase.set( 'Observer', Variable="CostFunctionJAtCurrentOptimum", String="print('\\n -----------------------------------\\n [VERBOSE] Current iteration: %i'%(len(var),))" )
513 __adaocase.set( 'Observer', Variable="CostFunctionJAtCurrentOptimum", Template="ValuePrinter", Info="\n [VERBOSE] Current best cost value:" )
514 __adaocase.set( 'Observer', Variable="CurrentOptimum", Template="ValuePrinter", Info="\n [VERBOSE] Current optimal state..:" )
518 InitialParameters = numpy.asarray(Background),
519 OptimalParameters = __adaocase.get("Analysis")[-1], # * Background,
520 OptimalSimulation = __adaocase.get("SimulatedObservationAtOptimum")[-1],
521 Measures = Observations,
522 NamesOfParameters = BNames,
523 NamesOfMeasures = ONames,
525 if self.__verbose: print("")
530 self.__result["SourceName"],
531 self.__result["Level"],
532 self.__result["Format"],
540 def calibrateMultiObs(self,
541 ModelName = None, # Model command file
542 ModelFormat = "Guess",
543 DataNames = None, # Multiple measures
544 DataFormat = "Guess",
545 LinkName = None, # CL with data names
546 LinkVariable = "Variable", # Name of variable names column
547 LinkFormat = "Guess",
548 BackgroundName = None, # Background
549 BackgroundFormat = "Guess",
550 MethodName = None, # Assimilation
551 MethodFormat = "Guess",
552 VariablesToCalibrate = None,
553 OutputVariables = None,
554 ResultName = "results.txt", # Results
555 ResultFormat = "Guess",
558 IntermediateInformation = False,
559 ComplexModel = False,
561 NeedInitVariables = False,
562 KeepCalculationFolders = False,
564 CreateOptimaldsin = False,
565 InitialSimulation = False,
566 VerifyFunction = False,
567 GlobalSensitivityCheck = False,
568 ParamSensitivityCheck = False,
569 CreateCSVoutputResults = False,
570 AdaptObservationError = False,
571 ResultsSummary = False,
572 AdvancedDebugModel = False,
573 TimeoutModelExecution = 300,
574 ModelStabilityCheck = False,
575 ModelStabilityCheckingPoints = 10
578 General method to set and realize a calibration (optimal estimation
579 task) of model parameters using multiple measures data, with an
580 explicit link between measures and data for the model.
583 self.__verbose = True
584 print(" [VERBOSE] Verbose mode activated")
585 if IntermediateInformation:
587 self.__IntermediateInformation = True
588 print(" [VERBOSE] IntermediateInformation provided")
590 self.__IntermediateInformation = False
591 print(" [VERBOSE] IntermediateInformation requires Verbose mode activated")
593 self.__IntermediateInformation = False
595 if TimeoutModelExecution < 0:
596 raise ValueError("TimeoutModelExecution must be positive")
598 if DataNames is not None:
599 # On force l'utilisateur à nommer dans son fichier de mesures
600 # les variables avec le même nom que des sorties dans le modèle
601 self.describeMultiMeasures(DataNames, OutputVariables, DataFormat)
604 print(" [VERBOSE] Timeoout for model execution is:", TimeoutModelExecution, "seconds")
606 if ModelFormat is not None: #included sot that the model format is allways in upper format
607 ModelFormat = ModelFormat.upper()
608 if ModelFormat == "DYMOLA": #to handle situations in which name Dymola is indicated
609 ModelFormat = "DYMOSIM"
610 if ModelName is not None:
611 self.describeModel(ModelName, VariablesToCalibrate, OutputVariables, ModelFormat)
612 if LinkName is not None:
613 # On force l'utilisateur à utiliser les fichiers de mesures requis
614 # et non pas tous ceux présents dans le Link
615 self.describeLink(DataNames, ModelName, LinkName, LinkVariable, LinkFormat)
616 if MethodName is not None:
617 self.describeMethod(MethodName, MethodFormat, BackgroundName, BackgroundFormat)
618 if ResultName is not None:
619 self.describeResult(ResultName, ResultLevel, ResultFormat)
621 #handling of the new option for complex models while keeping the previous key word in the code
623 print(" The option ComplexModel is deprecated and has to be replaced by NeedInitVariables. Please update your code.")
624 NeedInitVariables = ComplexModel
625 ComplexModel = NeedInitVariables
628 ONames, Observations = _readMultiMeasures(
629 self.__measures["SourceName"],
630 self.__measures["Variables"],
631 self.__measures["Format"],
635 Algo, Params, CovB, CovR = _readMethod(
636 self.__method["SourceName"],
637 self.__method["Format"],
640 if self.__method["BackgroundFormat"] == "DSIN": __bgfile = self.__model["SourceName"]
641 if self.__method["BackgroundFormat"] == "ADAO": __bgfile = self.__method["SourceName"]
642 if self.__method["BackgroundFormat"] == "USER": __bgfile = self.__method["BackgroundName"]
643 if self.__method["BackgroundFormat"] == "Guess":
644 if self.__method["BackgroundName"] is not None: __bgfile = self.__method["BackgroundName"]
645 if self.__method["SourceName"] is not None: __bgfile = self.__method["SourceName"]
646 if self.__model["SourceName"] is not None: __bgfile = self.__model["SourceName"]
647 BNames, Background, Bounds = _readBackground(
649 self.__model["VariablesToCalibrate"],
650 self.__method["BackgroundFormat"],
652 if "Bounds" not in Params: # On force la priorité aux bornes utilisateur
653 Params["Bounds"] = Bounds
654 LNames, LColumns, LVariablesToChange = _readLink(
655 self.__link["SourceName"],
656 self.__link["MeasuresNames"],
657 self.__link["Variable"],
658 self.__link["Format"],
661 print(" [VERBOSE] Measures read")
662 print(" [VERBOSE] Links read")
663 print(" [VERBOSE] Optimization information read")
664 print(" [VERBOSE] Background read:",Background)
665 print(" [VERBOSE] Bounds read:",Params['Bounds'])
667 if "DifferentialIncrement" not in Params:
668 Params["DifferentialIncrement"] = 0.001
669 if "EnableMultiProcessing" not in Params:
670 Params["EnableMultiProcessing"] = 0
671 if "NumberOfProcesses" not in Params:
672 Params["NumberOfProcesses"] = 0
673 if "StoreSupplementaryCalculations" not in Params:
674 Params["StoreSupplementaryCalculations"] = ["SimulatedObservationAtOptimum",]
675 if "StoreSupplementaryCalculations" in Params and \
676 "SimulatedObservationAtOptimum" not in Params["StoreSupplementaryCalculations"]:
677 Params["StoreSupplementaryCalculations"].append("SimulatedObservationAtOptimum")
678 if self.__verbose and "StoreSupplementaryCalculations" in Params and \
679 "CurrentOptimum" not in Params["StoreSupplementaryCalculations"]:
680 Params["StoreSupplementaryCalculations"].append("CurrentOptimum")
681 if self.__verbose and "StoreSupplementaryCalculations" in Params and \
682 "CostFunctionJAtCurrentOptimum" not in Params["StoreSupplementaryCalculations"]:
683 Params["StoreSupplementaryCalculations"].append("CostFunctionJAtCurrentOptimum")
685 dict_dsin_path_ref ={} #creation of the dict_dsin_ref (for storing the path of the dsfinal file for each dataset
688 VariablesToCalibrate = list(BNames)
690 def exedymosimMultiobs_REF_CALCULATION( x_values_matrix, dict_dsin_paths = dict_dsin_path_ref): #2ème argument à ne pas modifier
691 "Appel du modèle et restitution des résultats"
693 x_values = list(numpy.ravel(x_values_matrix))
694 x_inputs = dict(zip(VariablesToCalibrate,x_values))
695 dict_inputs_for_CWP = {}
697 # multi_y_values_matrix = []
698 for etat in range(len(LNames)):
699 simudir = self._setModelTmpDir_REF_CALCULATION(ModelName, ModelFormat, LNames[etat])
701 try: #to handle situations in which there is only one CL (boundary conditions) file
702 var_int = numpy.ravel(LColumns[:,etat])
704 var_int = numpy.ravel(LColumns)
705 dict_inputs = dict(zip(LVariablesToChange, var_int))
707 dict_inputs.update( x_inputs )
709 auto_simul = Automatic_Simulation(
711 dymosim_path = os.path.join(simudir, os.pardir, os.pardir) ,
712 dsres_path = os.path.join(simudir, str("dsres" +LNames[etat]+ ".mat")),
713 dict_inputs=dict_inputs,
715 timeout=TimeoutModelExecution,
716 without_modelicares=True,
718 auto_simul.single_simulation()
722 if AdvancedDebugModel:
723 dslog_file = os.path.join(simudir, 'dslog.txt')
724 debug_model_file = os.path.join(simudir, os.pardir, 'log_debug_model.txt')
726 shutil.copy(dslog_file,debug_model_file)
730 if auto_simul.success_code == 2 :
731 raise ValueError("The simulation falis after initialization, please check the log file and your model in order to be sure that the whole simulation can be performed (specially for dynamic models)")
733 if auto_simul.success_code == 3 :
734 raise ValueError("The simulation did not reach the final time, probably due to a lack of time, please increase the time for model execution by modifying the TimeoutModelExecution value in configuration.py")
736 if auto_simul.success_code == 0:
738 path_around_simu_no_CWP = os.path.join(simudir, os.pardir, os.pardir)
739 around_simu_no_CWP = Around_Simulation(dymola_version='2012',
740 curdir=path_around_simu_no_CWP,
741 source_list_iter_var = 'from_dymtranslog',
742 copy_from_dym_trans_log= os.path.join(path_around_simu_no_CWP,'ini.txt'),
743 mat = os.path.join(simudir, str("dsres" +LNames[etat]+ ".mat")),
744 iter_var_values_options={'source':'from_mat', 'moment':'initial'},
746 around_simu_no_CWP.set_list_iter_var()
747 around_simu_no_CWP.set_dict_iter_var()
749 writer_no_CWP = write_in_dsin.Write_in_dsin(dict_inputs =around_simu_no_CWP.dict_iter_var,filedir=simudir,dsin_name='dsin.txt',old_file_name='old_dsin.txt',new_file_name=str("dsin_updated_" +LNames[etat]+ ".txt"))
750 writer_no_CWP.write_in_dsin
751 writer_no_CWP.write_in_dsin()
753 dict_dsin_paths[LNames[etat]] = os.path.join(simudir,str("dsin_updated_" +LNames[etat]+ ".txt"))
756 from pst4mod.modelica_libraries.change_working_point import Dict_Var_To_Fix
757 from pst4mod.modelica_libraries.change_working_point import Working_Point_Modification
759 path_around_simu = os.path.join(simudir,os.path.pardir,os.path.pardir)
761 if not(os.path.exists(os.path.join(path_around_simu, str("ref" + ".mat")))):
762 auto_simul_ref = Automatic_Simulation(
763 simu_dir=path_around_simu,
764 dymosim_path = os.path.join(simudir, os.pardir, os.pardir) ,
767 dsres_path = os.path.join(path_around_simu, str("ref" + ".mat")),
768 timeout=TimeoutModelExecution,
769 without_modelicares=True,
771 auto_simul_ref.single_simulation()
773 if not(auto_simul_ref.success_code == 0):
774 raise ValueError("A working dsin.txt file must be provided, provide a dsin.txt that makes it possible for the initial simulation to converge")
776 temporary_files_to_remove = [os.path.join(path_around_simu,'dsfinal.txt'),
777 os.path.join(path_around_simu,'dsin_old.txt'),
778 os.path.join(path_around_simu,'dslog.txt'),
779 os.path.join(path_around_simu,'status'),
780 os.path.join(path_around_simu,'success')
782 for path_to_remove in temporary_files_to_remove:
783 if os.path.exists(path_to_remove):
784 os.remove(path_to_remove)
786 around_simu = Around_Simulation(dymola_version='2012',
787 curdir=path_around_simu,
788 source_list_iter_var = 'from_dymtranslog',
789 copy_from_dym_trans_log= os.path.join(path_around_simu,'ini.txt'),
790 mat = os.path.join(path_around_simu, str("ref" + ".mat")),
791 iter_var_values_options={'source':'from_mat'},
793 around_simu.set_list_iter_var()
794 around_simu.set_dict_iter_var()
796 reader_CWP = Reader(os.path.join(path_around_simu, str("ref" + ".mat")),'dymola')
797 x_values_intemediary = [reader_CWP.values(x_name)[1][-1] for x_name in VariablesToCalibrate]
798 x_values_from_ref_calculation = numpy.asarray(x_values_intemediary)
799 x_inputs_from_ref_calculation = dict(zip(VariablesToCalibrate,x_values_from_ref_calculation))
801 cl_values_intemediary = [reader_CWP.values(cl_name)[1][-1] for cl_name in LVariablesToChange]
802 cl_values_from_ref_calculation = numpy.asarray(cl_values_intemediary)
803 cl_inputs_from_ref_calculation = dict(zip(LVariablesToChange,cl_values_from_ref_calculation))
806 for var_calib in x_inputs_from_ref_calculation.keys():
807 dict_inputs_for_CWP[var_calib] = (x_inputs_from_ref_calculation[var_calib], x_inputs[var_calib])
809 try: #to handle situations in which there is only one CL (boundary conditions) file
810 var_int = numpy.ravel(LColumns[:,etat])
812 var_int = numpy.ravel(LColumns)
814 dict_inputs = dict(zip(LVariablesToChange, var_int))
816 for var_cl in cl_inputs_from_ref_calculation.keys():
817 dict_inputs_for_CWP[var_cl] = (cl_inputs_from_ref_calculation[var_cl], dict_inputs[var_cl])
819 dict_var_to_fix2 = Dict_Var_To_Fix(option='automatic',
820 dict_auto_var_to_fix=dict_inputs_for_CWP)
822 dict_var_to_fix2.set_dict_var_to_fix()
825 LOG_FILENAME = os.path.join(simudir,'change_working_point.log')
826 for handler in logging.root.handlers[:]:
827 logging.root.removeHandler(handler)
828 logging.basicConfig(filename=LOG_FILENAME,level=logging.DEBUG,filemode='w')
830 working_point_modif = Working_Point_Modification(main_dir = simudir,
831 simu_material_dir = simudir,
832 dymosim_path = os.path.join(simudir, os.pardir, os.pardir),
834 store_res_dir = 'RES',
835 dict_var_to_fix = dict_var_to_fix2.dict_var_to_fix,
836 around_simulation0 = around_simu,
837 var_to_follow_path= os.path.join(simudir,'var_to_follow.csv'),
838 gen_scripts_ini= False,
839 nit_max= 1000000000000000,
840 min_step_val = 0.00000000000005,
841 timeout = TimeoutModelExecution,
844 working_point_modif.create_working_directory()
845 working_point_modif.working_point_modification(skip_reference_simulation=True)
847 if AdvancedDebugModel:
848 dslog_file = os.path.join(simudir, 'SIMUS', 'dslog.txt')
849 debug_model_file = os.path.join(simudir, os.pardir, 'log_debug_model.txt')
851 shutil.copy(dslog_file,debug_model_file)
855 if not(os.path.exists(os.path.join(simudir, 'RES','1.0.mat'))):
856 raise ValueError(str("Simulation with Background values does not converge automatically, try to give a dsin.txt file that makes it possible to run a simulation with boundary conditions of " + LNames[etat]))
859 os.path.join(simudir, 'RES','1.0.mat'),
860 os.path.join(simudir, str("dsres" +LNames[etat]+ ".mat"))
863 path_around_simu_after_CWP = os.path.join(simudir,os.path.pardir,os.path.pardir)
864 around_simu_after_CWP = Around_Simulation(dymola_version='2012',
865 curdir=path_around_simu_after_CWP,
866 source_list_iter_var = 'from_dymtranslog',
867 copy_from_dym_trans_log= os.path.join(path_around_simu_after_CWP,'ini.txt'),
868 mat = os.path.join(simudir, 'RES','1.0.mat'), #Le fichier .mat qui vient d'être créé et pour lequel on atteint les valeurs cibles de la variable
869 iter_var_values_options={'source':'from_mat', 'moment':'initial'},
871 around_simu_after_CWP.set_list_iter_var()
872 around_simu_after_CWP.set_dict_iter_var()
874 writer_after_CWP = write_in_dsin.Write_in_dsin(dict_inputs =around_simu_after_CWP.dict_iter_var,filedir=os.path.join(simudir, 'SIMUS'),dsin_name='dsin.txt',old_file_name='old_dsin.txt',new_file_name='dsin_updated.txt')
875 writer_after_CWP.write_in_dsin
876 writer_after_CWP.write_in_dsin()
879 os.path.join(simudir, 'SIMUS','dsin_updated.txt'),
880 os.path.join(simudir, str("dsin_updated_" +LNames[etat]+ ".txt"))
882 dict_dsin_paths[LNames[etat]] = os.path.join(simudir,str("dsin_updated_" +LNames[etat]+ ".txt"))
885 os.path.join(os.path.dirname(dict_dsin_paths[LNames[etat]]), str("dsres" +LNames[etat]+ ".mat")),
886 os.path.join(os.path.dirname(dict_dsin_paths[LNames[etat]]), str("REF_for_CWP" + ".mat"))
893 def preparation_exedymosimMultiobs_simple():
894 from pst4mod.modelica_libraries import write_in_dsin
896 for etat in range(len(LNames)):
897 simudir = self._setModelTmpDir_simple(ModelName, ModelFormat, str("dsin_" +LNames[etat]+ ".txt"))
899 try: #to handle situations in which there is only one CL (boundary conditions) file
900 var_int = numpy.ravel(LColumns[:,etat])
902 var_int = numpy.ravel(LColumns)
903 dict_inputs = dict(zip(LVariablesToChange, var_int))
904 writer_no_CWP = write_in_dsin.Write_in_dsin(dict_inputs = dict_inputs, filedir=simudir, dsin_name=str("dsin_" +LNames[etat]+ ".txt"), old_file_name='old_dsin.txt',new_file_name = str("dsin_" +LNames[etat]+ ".txt"))
905 writer_no_CWP.write_in_dsin
906 writer_no_CWP.write_in_dsin()
909 def preparation_exefmu_om_Multiobs():
910 self._setModelTmpDir_simple(ModelName, ModelFormat, ModelName)
912 def exepythonMultiobs( x_values_matrix ):
913 "Appel du modèle et restitution des résultats"
915 x_values = list(numpy.ravel(x_values_matrix))
916 x_inputs = dict(zip(VariablesToCalibrate,x_values))
918 multi_y_values_matrix = []
919 for etat in range(len(LNames)):
920 simudir = self._setModelTmpDir(ModelName, ModelFormat)
922 dict_inputs = dict(zip(LVariablesToChange, numpy.ravel(LColumns[:,etat])))
923 dict_inputs.update( x_inputs )
925 auto_simul = Python_Simulation(
927 val_inputs=x_values_matrix,
929 timeout=TimeoutModelExecution,
930 verbose=self.__verbose)
931 y_values_matrix = auto_simul.single_simulation()
932 multi_y_values_matrix.append( y_values_matrix )
934 shutil.rmtree(simudir, ignore_errors=True)
935 y_values_matrix = numpy.ravel(numpy.array(multi_y_values_matrix))
936 return y_values_matrix
940 #Test for output variables - not repetition - This test must be before the redefinition of OutputVariables
941 OutputVariables_set = set(OutputVariables)
942 if len(OutputVariables) != len(OutputVariables_set):
943 diff_number_values = len(OutputVariables) - len(OutputVariables_set)
944 raise ValueError("There are " + str(diff_number_values) + " repeated output variables in the definition of OutputVariables, please remove repeated names.")
946 #Handle several times same obs -START
947 OutputVariables_new_tmp = []
949 whole_obs_in_fileobs = readObsnamesfile(self.__measures["SourceName"])
951 for outputname in OutputVariables:
952 if outputname not in whole_obs_in_fileobs:
953 print(" [VERBOSE] /!\\ ", outputname," is not defined in the measurements file")
954 for obsname in whole_obs_in_fileobs:
955 if obsname in OutputVariables:
956 OutputVariables_new_tmp.append(obsname)
958 OutputVariables_new = []
959 for obs_name in OutputVariables: #to get the same order as for measures
960 count_obs_name = OutputVariables_new_tmp.count(obs_name)
961 for i in range(count_obs_name):
962 OutputVariables_new.append(obs_name)
963 OutputVariables = OutputVariables_new
965 ONames = OutputVariables
966 #Handle several times same obs - END
969 ODates,OHours = readMultiDatesHours(self.__measures["SourceName"])
971 List_Multideltatime =[]
973 for etat in range(len(ODates)):
974 time_initial = datetime.strptime(ODates[etat][0]+ ' ' + OHours[etat][0], '%d/%m/%Y %H:%M:%S')
977 if len(ODates[etat])>1:
978 list_deltatime.append(0)
979 for index_tmp in range(len(ODates[etat])-1): #the first date is considered as reference
980 time_index_tmp = datetime.strptime(ODates[etat][index_tmp+1] + ' ' + OHours[etat][index_tmp+1], '%d/%m/%Y %H:%M:%S')
981 list_deltatime.append((time_index_tmp-time_initial).total_seconds())
983 if (time_index_tmp-time_initial).total_seconds() < 0:
984 raise ValueError('The initial date is not chronological for state %s'%str(etat))
986 List_Multideltatime.append(list_deltatime)
991 number_observations_by_cl = [] #to be used in the test function, number of measurements by boundary condition (i.e. number of lines in the measurements file)
993 for etat in range(len(ODates)): #correspond au nombre de conditions aux limites
994 Obs_etat = Observations[etat].tolist()
996 for obs in range(len(list(ONames))):
997 if len(ODates[etat])==1: #cas classique, statique
998 list_obs_etat = Obs_etat
999 if not(isinstance(list_obs_etat, list)):
1000 list_obs_etat = [list_obs_etat]
1003 for time_etat in range(len(Obs_etat)):
1004 if not(isinstance(Obs_etat[time_etat], list)): #une seule grandeur observée
1005 list_obs_etat.append(Obs_etat[time_etat])
1007 list_obs_etat.append(Obs_etat[time_etat][obs])
1009 #EXTRA for verify tests - START
1010 if len(ODates[etat])==1:
1011 number_observations_by_cl.append(1) #in this case only one line in the measurements file, this is just to be used as input of the test functions
1013 number_observations_by_cl.append(len(Obs_etat)) #we get the number of lines of the measurements file, this is just to be used as input of the test functions
1014 #EXTRA for verify tests - END
1016 Observations_new = Observations_new + list_obs_etat
1017 Observations = Observations_new
1023 def prev_adao_version_TOP_LEVEL_exedymosimMultiobs( x_values_matrix):
1024 y_values_matrix = TOP_LEVEL_exedymosimMultiobs( x_values_matrix=x_values_matrix, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, ModelFormat = ModelFormat, KeepCalculationFolders = KeepCalculationFolders, Verbose = Verbose, dict_dsin_paths = dict_dsin_path_ref, Linux = Linux, List_Multideltatime = List_Multideltatime, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution)
1025 return y_values_matrix
1027 def prev_adao_version_TOPLEVEL_exedymosimMultiobs_simple(x_values_matrix):
1028 y_values_matrix = TOPLEVEL_exedymosimMultiobs_simple(x_values_matrix = x_values_matrix, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), KeepCalculationFolders = KeepCalculationFolders, Linux = Linux, List_Multideltatime = List_Multideltatime, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution)
1029 return y_values_matrix
1031 if self.__model["Format"].upper() in ["DYMOSIM", "GUESS"]:
1033 if adao.version[:5] < '9.7.0' and adao.version[5]=='.': #to handle differences in the way directoperator is mangaged after modifications brought to direct operator
1035 simulateur = prev_adao_version_TOP_LEVEL_exedymosimMultiobs
1036 exedymosimMultiobs_REF_CALCULATION(Background)
1038 simulateur = prev_adao_version_TOPLEVEL_exedymosimMultiobs_simple
1039 preparation_exedymosimMultiobs_simple()
1043 simulateur = TOP_LEVEL_exedymosimMultiobs #to handle parallel computing
1044 exedymosimMultiobs_REF_CALCULATION(Background)
1046 simulateur = TOPLEVEL_exedymosimMultiobs_simple #to handle parallel computing
1047 preparation_exedymosimMultiobs_simple()
1049 elif self.__model["Format"].upper() == "PYSIM":
1050 simulateur = exepythonMultiobs
1051 elif self.__model["Format"].upper() == "FMI" or self.__model["Format"].upper() == "FMU":
1052 simulateur = TOP_LEVEL_exefmuMultiobs
1053 preparation_exefmu_om_Multiobs()
1054 elif self.__model["Format"].upper() == "OPENMODELICA":
1055 preparation_exefmu_om_Multiobs()
1056 simulateur = TOP_LEVEL_exeOpenModelicaMultiobs
1059 raise ValueError("Unknown model format \"%s\""%self.__model["Format"].upper())
1061 def verify_function(x_values, model_format):
1064 simulation_results_all = []
1065 list_simulation_time = []
1066 for i in list(range(3)): #we consider only three iterations
1067 start_time_verify = time.time()
1068 if model_format in ["DYMOSIM"]:
1070 simulation_results = prev_adao_version_TOP_LEVEL_exedymosimMultiobs(x_values)
1072 simulation_results = prev_adao_version_TOPLEVEL_exedymosimMultiobs_simple(x_values)
1073 elif model_format in ["FMI","FMU"]:
1074 simulation_results = TOP_LEVEL_exefmuMultiobs(x_values, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution, FMUInput = FMUInput)
1075 elif model_format in ["OPENMODELICA"]:
1076 simulation_results = TOP_LEVEL_exeOpenModelicaMultiobs(x_values, KeepCalculationFolders = KeepCalculationFolders, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, Linux = Linux, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution)
1078 raise NotImplementedError("Not yet implemented for current model format: ", model_format)
1079 end_time_verify = time.time()
1081 list_simulation_time.append(round(end_time_verify - start_time_verify,3))
1082 simulation_results_all.append(simulation_results)
1085 elapsed_time = round(sum(list_simulation_time)/len(list_simulation_time),3)
1086 #not necessary, already in good format
1087 # numpy.set_printoptions(precision=2,linewidth=5000, threshold = 10000)
1088 # reshaped_simulation_results_all = numpy.array(simulation_results_all).reshape((len(OutputVariables),-1)).transpose()
1090 # numpy.set_printoptions(precision=2,linewidth=5000)
1092 list_check_output = []
1093 list_false_output = []
1094 simulation_results_all_one_boundary_condition =[]
1096 for simulation_output in simulation_results_all: #iteration over the three simulations performed
1097 # THIS DOES NOT WORK IF OBSERVATIONS FILES ARE OF DIFFERENT SIZE# simultation_ouput_reshaped = simulation_output.reshape((len(LNames),-1))[0] #consideration of only one set of boundary conditions (otherwise outputs could be modified by boundary condition modifications)
1098 total_number_measurements_first_simulation = number_observations_by_cl[0]*len(OutputVariables)
1099 first_simulation_results = simulation_output[:total_number_measurements_first_simulation]
1100 simulation_results_all_one_boundary_condition.append(first_simulation_results)
1102 for index_line_time_simulation in range(len(simulation_results_all_one_boundary_condition[0].reshape(len(OutputVariables),-1).transpose())): #iteration over each time results individually
1103 for i in range(len(OutputVariables)):
1104 if (simulation_results_all_one_boundary_condition[0].reshape(len(OutputVariables),-1).transpose()[index_line_time_simulation][i] == simulation_results_all_one_boundary_condition[1].reshape(len(OutputVariables),-1).transpose()[index_line_time_simulation][i]) & (simulation_results_all_one_boundary_condition[0].reshape(len(OutputVariables),-1).transpose()[index_line_time_simulation][i] == simulation_results_all_one_boundary_condition[2].reshape(len(OutputVariables),-1).transpose()[index_line_time_simulation][i]):
1105 list_check_output.append("True")
1108 list_check_output.append("False")
1109 if OutputVariables[i] in list_false_output: #problematic variable already included
1112 list_false_output.append(OutputVariables[i])
1114 if "False" in list_check_output:
1115 verify_function_result = "False"
1117 verify_function_result = "True"
1118 return verify_function_result, list_false_output, elapsed_time
1120 def verify_function_sensitivity(x_values, increment, model_format):
1122 if model_format in ["DYMOSIM"]:
1124 simulation_results = prev_adao_version_TOP_LEVEL_exedymosimMultiobs(x_values)
1126 simulation_results = prev_adao_version_TOPLEVEL_exedymosimMultiobs_simple(x_values)
1127 elif model_format in ["FMI","FMU"]:
1128 simulation_results = TOP_LEVEL_exefmuMultiobs(x_values, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution, FMUInput = FMUInput)
1129 elif model_format in ["OPENMODELICA"]:
1130 simulation_results = TOP_LEVEL_exeOpenModelicaMultiobs(x_values, KeepCalculationFolders = KeepCalculationFolders, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, Linux = Linux, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution)
1132 raise NotImplementedError("Not yet implemented for current model format: ", model_format)
1135 x_values_modif = [x*(1+increment) for x in x_values]
1137 if model_format in ["DYMOSIM"]:
1139 simulation_results_modif = prev_adao_version_TOP_LEVEL_exedymosimMultiobs(x_values_modif)
1141 simulation_results_modif = prev_adao_version_TOPLEVEL_exedymosimMultiobs_simple(x_values_modif)
1142 elif model_format in ["FMI","FMU"]:
1143 simulation_results_modif = TOP_LEVEL_exefmuMultiobs(x_values_modif, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution, FMUInput = FMUInput)
1144 elif model_format in ["OPENMODELICA"]:
1145 simulation_results_modif = TOP_LEVEL_exeOpenModelicaMultiobs(x_values_modif, KeepCalculationFolders = KeepCalculationFolders, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, Linux = Linux, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution)
1147 raise NotImplementedError("Not yet implemented for current model format: ", model_format)
1149 list_check_output = ["False" for i in OutputVariables] #everything initialized as not good and then if something changes (at least one point), then it is OK
1150 list_not_modified_variables =[]
1152 total_number_measurements_first_simulation = number_observations_by_cl[0]*len(OutputVariables)
1154 # THIS DOES NOT WORK IF OBSERVATIONS FILES ARE OF DIFFERENT SIZE# simulation_results_one_boundary_condition = simulation_results.reshape((len(LNames),-1))[0]
1155 # THIS DOES NOT WORK IF OBSERVATIONS FILES ARE OF DIFFERENT SIZE# simulation_results_modif_one_boundary_condition = simulation_results_modif.reshape((len(LNames),-1))[0]
1156 simulation_results_one_boundary_condition = simulation_results[:total_number_measurements_first_simulation]
1157 simulation_results_modif_one_boundary_condition = simulation_results_modif[:total_number_measurements_first_simulation]
1159 for index_line_time_simulation in range(len(simulation_results_one_boundary_condition.reshape(len(OutputVariables),-1).transpose())):
1161 for i in range(len(OutputVariables)):
1162 if simulation_results_one_boundary_condition.reshape(len(OutputVariables),-1).transpose()[index_line_time_simulation][i] != simulation_results_modif_one_boundary_condition.reshape(len(OutputVariables),-1).transpose()[index_line_time_simulation][i] :
1163 list_check_output[i] = "True"
1166 for i in range(len(OutputVariables)):
1167 if list_check_output[i] == "False": #it means the variable has not changed
1168 list_not_modified_variables.append(OutputVariables[i])
1170 if "False" in list_check_output:
1171 verify_function_sensitivity_result = "False"
1173 verify_function_sensitivity_result = "True"
1175 return verify_function_sensitivity_result, list_not_modified_variables
1177 def verify_function_sensitivity_parameter(x_values, increment, model_format):
1179 if model_format in ["DYMOSIM"]:
1181 simulation_results_ref = prev_adao_version_TOP_LEVEL_exedymosimMultiobs(x_values)
1183 simulation_results_ref = prev_adao_version_TOPLEVEL_exedymosimMultiobs_simple(x_values)
1184 elif model_format in ["FMI","FMU"]:
1185 simulation_results_ref = TOP_LEVEL_exefmuMultiobs(x_values, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution, FMUInput = FMUInput)
1186 elif model_format in ["OPENMODELICA"]:
1187 simulation_results_ref = TOP_LEVEL_exeOpenModelicaMultiobs(x_values, KeepCalculationFolders = KeepCalculationFolders, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, Linux = Linux, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution)
1189 raise NotImplementedError("Not yet implemented for current model format: ", model_format)
1191 list_param_with_no_impact = []
1192 Check_params_impact = "True"
1194 x_values_modif = [x*(1+increment) for x in x_values]
1197 for param_index in range(len(x_values)): #the goal is to modify each param individually and check if it is possible to observe variations
1198 x_values_param = list(x_values) #avoid memory surprises
1199 x_values_param[param_index] = x_values_modif[param_index]
1200 if model_format in ["DYMOSIM"]:
1202 simulation_results_param = prev_adao_version_TOP_LEVEL_exedymosimMultiobs(x_values_param)
1204 simulation_results_param = prev_adao_version_TOPLEVEL_exedymosimMultiobs_simple(x_values_param)
1205 elif model_format in ["FMI","FMU"]:
1206 simulation_results_param = TOP_LEVEL_exefmuMultiobs(x_values_param, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution, FMUInput = FMUInput)
1207 elif model_format in ["OPENMODELICA"]:
1208 simulation_results_param = TOP_LEVEL_exeOpenModelicaMultiobs(x_values_param, KeepCalculationFolders = KeepCalculationFolders, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, Linux = Linux, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution)
1210 raise NotImplementedError("Not yet implemented for current model format: ", model_format)
1212 if numpy.array_equal(simulation_results_param,simulation_results_ref): #there is not impact on the whole simulation
1213 list_param_with_no_impact.append(VariablesToCalibrate[param_index])
1214 Check_params_impact = "False"
1216 return Check_params_impact, list_param_with_no_impact
1219 def check_model_stability(x_values_names, background_values, bounds, number_of_tests, model_format):
1221 dict_check_model_stability = {}
1222 for index_x in range(len(x_values_names)):
1223 x_bg_tested = x_values_names[index_x]
1225 print(" [VERBOSE] Model stability check is being performed for the following variable to be calibrated: ", x_bg_tested)
1227 if bounds[index_x][0] == None or bounds[index_x][0] == "None":
1228 bounds_inf_x_bg_tested = background_values[index_x] - abs(background_values[index_x]*100)
1230 bounds_inf_x_bg_tested = bounds[index_x][0]
1232 if bounds[index_x][1] == None or bounds[index_x][1] == "None":
1233 bounds_sup_x_bg_tested = background_values[index_x] + abs(background_values[index_x]*100)
1235 bounds_sup_x_bg_tested = bounds[index_x][1]
1237 #avoid problems with values to be tested hen the bounds are strictly equal to 0
1238 if bounds_inf_x_bg_tested == 0: #which means that bounds_sup_x_bg_tested > 0
1239 bounds_inf_x_bg_tested = min(1e-9,bounds_sup_x_bg_tested*1e-9)
1240 if bounds_sup_x_bg_tested == 0: #which means that bounds_inf_x_bg_tested < 0
1241 bounds_sup_x_bg_tested = max(-1e-9,bounds_inf_x_bg_tested*1e-9)
1243 list_x_bg_tested = []
1245 if bounds_inf_x_bg_tested < 0:
1246 if bounds_sup_x_bg_tested<0:
1247 list_x_bg_tested = list(numpy.geomspace(bounds_inf_x_bg_tested,bounds_sup_x_bg_tested,int(number_of_tests)))
1249 list_x_bg_tested = list(numpy.geomspace(bounds_inf_x_bg_tested,bounds_inf_x_bg_tested*1e-4,int(number_of_tests/2)))
1250 list_x_bg_tested = list_x_bg_tested + list((numpy.geomspace(bounds_sup_x_bg_tested*1e-4,bounds_sup_x_bg_tested,int(number_of_tests/2))))
1251 else: #both bounds are >0
1252 list_x_bg_tested = list(numpy.geomspace(bounds_inf_x_bg_tested,bounds_sup_x_bg_tested,int(number_of_tests)))
1254 #for a linear partition, we consider a log repartition better
1255 # for index_partition in range(number_of_tests):
1256 # list_x_bg_tested.append(bounds_inf_x_bg_tested + index_partition*(bounds_sup_x_bg_tested-bounds_inf_x_bg_tested)/(number_of_tests-1))
1258 list_failed_model_evaluation = []
1259 for x_value_bg_tested in list_x_bg_tested:
1261 x_whole_values_with_bg = copy.deepcopy(background_values) #avoid modifiction of backgroundvalues
1262 x_whole_values_with_bg[index_x] = x_value_bg_tested
1264 if model_format in ["DYMOSIM"]:
1266 simulation_results = prev_adao_version_TOP_LEVEL_exedymosimMultiobs(x_whole_values_with_bg)
1268 simulation_results = prev_adao_version_TOPLEVEL_exedymosimMultiobs_simple(x_whole_values_with_bg)
1269 elif model_format in ["FMI","FMU"]:
1270 simulation_results = TOP_LEVEL_exefmuMultiobs(x_whole_values_with_bg, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution, FMUInput = FMUInput)
1271 elif model_format in ["OPENMODELICA"]:
1272 simulation_results = TOP_LEVEL_exeOpenModelicaMultiobs(x_whole_values_with_bg, KeepCalculationFolders = KeepCalculationFolders, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, Linux = Linux, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution)
1274 raise NotImplementedError("Not yet implemented for current model format: ", model_format)
1276 list_failed_model_evaluation.append(x_value_bg_tested)
1278 dict_check_model_stability[x_bg_tested] = list_failed_model_evaluation
1279 return dict_check_model_stability
1282 check_function, list_not_same_model_outputs, elapsed_time_for_simu = verify_function(x_values = Background, model_format = self.__model["Format"].upper())
1283 if check_function == "True":
1284 print(" [VERBOSE] Results of function checking: ", "OK for all model outputs, ","simulation time (including all boundary conditions) is ",elapsed_time_for_simu, " seconds" )
1286 print(" [VERBOSE] Results of function checking: ", "NOT OK for all model outputs")
1287 print(" [VERBOSE] The following model outputs are not the same when repeating a simulation: ", list_not_same_model_outputs )
1289 if GlobalSensitivityCheck:
1290 check_sensitivity, list_same_model_outputs = verify_function_sensitivity(x_values = Background, increment = Params["DifferentialIncrement"], model_format = self.__model["Format"].upper())
1291 if VerifyFunction == False:
1292 print(" [VERBOSE] /!\\ Activate VerifyFunction option and check that the result of the check is OK to get reliable results from GlobalSensitivityCheck option")
1293 if check_sensitivity == "True":
1294 print(" [VERBOSE] Results of function global sensitivity analysis: ", "All model outputs vary when background values are modified")
1296 print(" [VERBOSE] Results of function global sensitivity analysis: ", "Some model outputs DO NOT vary when background values are modified")
1297 print(" [VERBOSE] The following model outputs DO NOT vary when background values are modified: ", list_same_model_outputs )
1299 if ParamSensitivityCheck:
1300 check_function, list_params_without_impact = verify_function_sensitivity_parameter(x_values = Background, increment = Params["DifferentialIncrement"], model_format = self.__model["Format"].upper())
1301 if VerifyFunction == False:
1302 print(" [VERBOSE] /!\\ Activate VerifyFunction option and check that the result of the check is OK to get reliable results from ParamSensitivityCheck option")
1303 if check_function == "True":
1304 print(" [VERBOSE] Results of parameter sensitivity checking: ", "OK, all parameters have an impact on model outputs")
1306 print(" [VERBOSE] Results of parameter sensitivity checking: ", "NOT OK, some parameters do not have an impact on model outputs")
1307 print(" [VERBOSE] The following parameters do not have an impact on model outputs: ", list_params_without_impact )
1309 if ModelStabilityCheck:
1310 print(" [VERBOSE] Model stability check is being performed, this check might take several minutes depending on the model and number of variables to calibrate (you can reduce the value of ModelStabilityCheckingPoints option to reduce checking time, by default its value is 10)")
1311 check_stability = check_model_stability(x_values_names = VariablesToCalibrate, background_values = Background, bounds = Params['Bounds'], number_of_tests = ModelStabilityCheckingPoints, model_format = self.__model["Format"].upper())
1312 if check_stability == dict(zip(VariablesToCalibrate,[list() for i in VariablesToCalibrate])):
1313 print(" [VERBOSE] Results of Model stability check: ", "The model simulates correctly when modifying individually the values of the variables to calibrate within the range in which their optimal value should be found")
1315 print(" [VERBOSE] Results of Model stability check: ", "The model FAILS to converge with the following values of variables to calibrate", check_stability)
1316 print(" [VERBOSE] The data assimilation procedure might not succeed if the model fails to converge, please revise the model and/or the range in which the optimal value of the variables to calibrate should be found")
1320 __adaocase = adaoBuilder.New()
1321 __adaocase.set( 'AlgorithmParameters', Algorithm=Algo, Parameters=Params )
1322 __adaocase.set( 'Background', Vector=Background)
1323 __adaocase.set( 'Observation', Vector=Observations )
1324 if type(CovB) is float:
1325 __adaocase.set( 'BackgroundError', ScalarSparseMatrix=CovB )
1327 __adaocase.set( 'BackgroundError', Matrix=CovB )
1329 if AdaptObservationError:
1330 Square_Observations = [x*x for x in Observations]
1332 # if (type(CovR) is float and CovR!= 1.0 ):# by default we still consider the square, better not to modify it
1333 if type(CovR) is float:
1334 ObsError = CovR*numpy.diag(Square_Observations)
1335 __adaocase.set( 'ObservationError', Matrix = ObsError )
1337 ObsError = 1e-15*numpy.diag(Square_Observations) #arbitary low value to eventually handle pressures in Pa
1338 __adaocase.set( 'ObservationError', Matrix = ObsError )
1339 print(" [VERBOSE] AdaptObservationError was set to True, the following ObservationError matrix has been considered: ")
1341 else: #classical situation
1342 if type(CovR) is float:
1343 __adaocase.set( 'ObservationError', ScalarSparseMatrix=CovR )
1345 __adaocase.set( 'ObservationError', Matrix=CovR )
1347 if Params["EnableMultiProcessing"]==1 :
1348 ParallelComputationGradient = True
1350 ParallelComputationGradient = False
1353 if self.__model["Format"].upper() in ["DYMOSIM"]:
1354 if ParallelComputationGradient:
1356 if adao.version[:5] < '9.7.0' and adao.version[5]=='.':
1357 raise ValueError("ADAO version must be at least 9.7.0 to support parallel computation of the gradient, current version is %s"%str(adao.version[:5]))
1360 __adaocase.set( 'ObservationOperator', OneFunction=simulateur,
1361 Parameters = {"DifferentialIncrement":Params["DifferentialIncrement"], "EnableMultiProcessing":Params["EnableMultiProcessing"], "NumberOfProcesses" : Params["NumberOfProcesses"]},
1362 ExtraArguments = {'VariablesToCalibrate':VariablesToCalibrate, 'OutputVariables' : OutputVariables, 'LNames': LNames, 'ModelFormat': ModelFormat, 'KeepCalculationFolders': KeepCalculationFolders, 'Verbose' : Verbose, 'dict_dsin_paths' : dict_dsin_path_ref, 'Linux': Linux, 'List_Multideltatime' : List_Multideltatime, 'AdvancedDebugModel':AdvancedDebugModel, 'TimeoutModelExecution':TimeoutModelExecution },
1365 __adaocase.set( 'ObservationOperator', OneFunction=simulateur,
1366 Parameters = {"DifferentialIncrement":Params["DifferentialIncrement"], "EnableMultiProcessing":Params["EnableMultiProcessing"], "NumberOfProcesses" : Params["NumberOfProcesses"]},
1367 ExtraArguments = {'VariablesToCalibrate':VariablesToCalibrate, 'OutputVariables' : OutputVariables, 'LNames': LNames, 'ref_simudir': self._get_Name_ModelTmpDir_simple(ModelName), 'KeepCalculationFolders': KeepCalculationFolders, 'Linux': Linux, 'List_Multideltatime' : List_Multideltatime, 'AdvancedDebugModel':AdvancedDebugModel, 'TimeoutModelExecution':TimeoutModelExecution },
1371 if adao.version[:5] < '9.7.0' and adao.version[5]=='.':
1372 __adaocase.set( 'ObservationOperator', OneFunction=simulateur, Parameters = {"DifferentialIncrement":Params["DifferentialIncrement"]})
1376 __adaocase.set( 'ObservationOperator', OneFunction=simulateur,
1377 Parameters = {"DifferentialIncrement":Params["DifferentialIncrement"]},
1378 ExtraArguments = {'VariablesToCalibrate':VariablesToCalibrate, 'OutputVariables' : OutputVariables, 'LNames': LNames, 'ModelFormat': ModelFormat, 'KeepCalculationFolders': KeepCalculationFolders, 'Verbose' : Verbose, 'dict_dsin_paths' : dict_dsin_path_ref, 'Linux': Linux, 'List_Multideltatime' : List_Multideltatime, 'AdvancedDebugModel':AdvancedDebugModel, 'TimeoutModelExecution':TimeoutModelExecution },
1381 __adaocase.set( 'ObservationOperator', OneFunction=simulateur,
1382 Parameters = {"DifferentialIncrement":Params["DifferentialIncrement"]},
1383 ExtraArguments = {'VariablesToCalibrate':VariablesToCalibrate, 'OutputVariables' : OutputVariables, 'LNames': LNames, 'ref_simudir': self._get_Name_ModelTmpDir_simple(ModelName), 'KeepCalculationFolders': KeepCalculationFolders, 'Linux': Linux, 'List_Multideltatime' : List_Multideltatime, 'AdvancedDebugModel':AdvancedDebugModel, 'TimeoutModelExecution':TimeoutModelExecution })
1385 elif self.__model["Format"].upper() in ["OPENMODELICA"]: #OPENMODELICA (different from FMI since there are different keywords for the function: Linux and KeepCalculationFolders)
1388 print("ComplexModel option is only valid for DYMOSIM format (it has no effect for the current calculation)")
1390 if ParallelComputationGradient:
1391 if adao.version[:5] < '9.7.0' and adao.version[5]=='.':
1392 raise ValueError("ADAO version must be at least 9.7.0 to support parallel computation of the gradient, current version is %s"%str(adao.version[:5]))
1394 __adaocase.set( 'ObservationOperator', OneFunction=simulateur,
1395 Parameters = {"DifferentialIncrement":Params["DifferentialIncrement"]},
1396 ExtraArguments = {'VariablesToCalibrate':VariablesToCalibrate, 'OutputVariables' : OutputVariables, 'LNames': LNames, 'LColumns':LColumns, 'LVariablesToChange':LVariablesToChange, 'ref_simudir': self._get_Name_ModelTmpDir_simple(ModelName), 'ModelName': ModelName, 'List_Multideltatime' : List_Multideltatime, 'Linux': Linux, 'KeepCalculationFolders': KeepCalculationFolders, 'AdvancedDebugModel':AdvancedDebugModel, 'TimeoutModelExecution':TimeoutModelExecution})
1399 if adao.version[:5] < '9.7.0' and adao.version[5]=='.':
1400 raise ValueError("ADAO version must be at least 9.7.0 to support other formats than DYMOSIM")
1402 __adaocase.set( 'ObservationOperator', OneFunction=simulateur,
1403 Parameters = {"DifferentialIncrement":Params["DifferentialIncrement"]},
1404 ExtraArguments = {'VariablesToCalibrate':VariablesToCalibrate, 'OutputVariables' : OutputVariables, 'LNames': LNames, 'LColumns':LColumns, 'LVariablesToChange':LVariablesToChange, 'ref_simudir': self._get_Name_ModelTmpDir_simple(ModelName), 'ModelName': ModelName, 'List_Multideltatime' : List_Multideltatime, 'Linux': Linux, 'KeepCalculationFolders': KeepCalculationFolders, 'AdvancedDebugModel':AdvancedDebugModel, 'TimeoutModelExecution':TimeoutModelExecution})
1406 else: #FMI or GUESS format
1409 print("ComplexModel option is only valid for DYMOSIM format (it has no effect for the current calculation)")
1411 if ParallelComputationGradient:
1412 if adao.version[:5] < '9.7.0' and adao.version[5]=='.':
1413 raise ValueError("ADAO version must be at least 9.7.0 to support parallel computation of the gradient, current version is %s"%str(adao.version[:5]))
1415 __adaocase.set( 'ObservationOperator', OneFunction=simulateur,
1416 Parameters = {"DifferentialIncrement":Params["DifferentialIncrement"]},
1417 ExtraArguments = {'VariablesToCalibrate':VariablesToCalibrate, 'OutputVariables' : OutputVariables, 'LNames': LNames, 'LColumns':LColumns, 'LVariablesToChange':LVariablesToChange, 'ref_simudir': self._get_Name_ModelTmpDir_simple(ModelName), 'ModelName': ModelName, 'List_Multideltatime' : List_Multideltatime, 'AdvancedDebugModel':AdvancedDebugModel, 'TimeoutModelExecution':TimeoutModelExecution, 'FMUInput':FMUInput})
1420 if adao.version[:5] < '9.7.0' and adao.version[5]=='.':
1421 raise ValueError("ADAO version must be at least 9.7.0 to support other formats than DYMOSIM")
1423 __adaocase.set( 'ObservationOperator', OneFunction=simulateur,
1424 Parameters = {"DifferentialIncrement":Params["DifferentialIncrement"]},
1425 ExtraArguments = {'VariablesToCalibrate':VariablesToCalibrate, 'OutputVariables' : OutputVariables, 'LNames': LNames, 'LColumns':LColumns, 'LVariablesToChange':LVariablesToChange, 'ref_simudir': self._get_Name_ModelTmpDir_simple(ModelName), 'ModelName': ModelName, 'List_Multideltatime' : List_Multideltatime, 'AdvancedDebugModel':AdvancedDebugModel, 'TimeoutModelExecution':TimeoutModelExecution, 'FMUInput':FMUInput})
1430 __adaocase.set( 'Observer', Variable="CostFunctionJAtCurrentOptimum", String="print('\\n -----------------------------------\\n [VERBOSE] Current iteration number: %i'%len(var))" )
1431 __adaocase.set( 'Observer', Variable="CostFunctionJAtCurrentOptimum", Template="ValuePrinter", Info="\n [VERBOSE] Current optimal cost value.:" )
1432 __adaocase.set( 'Observer', Variable="CurrentOptimum", Template="ValuePrinter", Info="\n [VERBOSE] Current optimal state......:" ) #
1434 if self.__IntermediateInformation:
1435 __adaocase.set( 'Observer', Variable="CostFunctionJ", Template="ValuePrinter", Info="\n [VERBOSE] Current cost value.:" )
1436 __adaocase.set( 'Observer', Variable="CurrentState", Template="ValuePrinter", Info="\n [VERBOSE] Current state......:" )
1438 __adaocase.execute()
1441 if InitialSimulation: #the prev_adao_version_TOP_LEVEL_XXX is already defined with proper arguments (boundary conditions etc.)
1443 if self.__model["Format"].upper() in ["DYMOSIM"]:
1445 initialsimulation_results = prev_adao_version_TOP_LEVEL_exedymosimMultiobs(Background)
1447 initialsimulation_results = prev_adao_version_TOPLEVEL_exedymosimMultiobs_simple(Background)
1448 elif self.__model["Format"].upper() in ["FMI", "FMU"]:
1449 initialsimulation_results = TOP_LEVEL_exefmuMultiobs(Background, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution, FMUInput = FMUInput)
1451 elif self.__model["Format"].upper() in ["OPENMODELICA"]:
1452 initialsimulation_results = TOP_LEVEL_exeOpenModelicaMultiobs(Background, KeepCalculationFolders = KeepCalculationFolders, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, Linux = Linux, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution )
1455 InitialParameters = numpy.asarray(Background),
1456 OptimalParameters = __adaocase.get("Analysis")[-1], # * Background,
1457 InitialSimulation = numpy.asarray(initialsimulation_results), #the simulation results are given for backgroung values, can be used to evaluate the overall performance of the calibration procedure
1458 OptimalSimulation = __adaocase.get("SimulatedObservationAtOptimum")[-1],
1459 Measures = Observations,
1460 NamesOfParameters = BNames,
1461 NamesOfMeasures = ONames,
1465 InitialParameters = numpy.asarray(Background),
1466 OptimalParameters = __adaocase.get("Analysis")[-1], # * Background,
1467 OptimalSimulation = __adaocase.get("SimulatedObservationAtOptimum")[-1],
1468 Measures = Observations,
1469 NamesOfParameters = BNames,
1470 NamesOfMeasures = ONames,
1473 if "APosterioriCovariance" in Params["StoreSupplementaryCalculations"]:
1474 __resultats["APosterioriCovariance"] = __adaocase.get("APosterioriCovariance")[-1]
1477 if InitialSimulation:
1478 measures_input = __resultats ["Measures"]
1479 initial_simu = __resultats ["InitialSimulation"]
1480 optimal_simu = __resultats ["OptimalSimulation"]
1482 diff_measures_initial_rmse = numpy.array([(measures_input[i]-initial_simu[i])**2 for i in range(len(measures_input))])
1483 diff_measures_initial_reldiff = numpy.array([abs((measures_input[i]-initial_simu[i])/measures_input[i]) for i in range(len(measures_input)) if measures_input[i] !=0])
1485 diff_measures_optimal_rmse = numpy.array([(measures_input[i]-optimal_simu[i])**2 for i in range(len(measures_input))])
1486 diff_measures_optimal_reldiff = numpy.array([abs((measures_input[i]-optimal_simu[i])/measures_input[i]) for i in range(len(measures_input)) if measures_input[i] !=0])
1488 diff_measures_initial_rmse_res = str(numpy.format_float_scientific(numpy.sqrt(diff_measures_initial_rmse.mean()), precision = 3))
1489 diff_measures_initial_reldiff_res = str(numpy.round(diff_measures_initial_reldiff.mean()*100,decimals =2))
1491 diff_measures_optimal_rmse_res = str(numpy.format_float_scientific(numpy.sqrt(diff_measures_optimal_rmse.mean()), precision = 3))
1492 diff_measures_optimal_reldiff_res = str(numpy.round(diff_measures_optimal_reldiff.mean()*100,decimals =2))
1494 __resultats["ResultsSummary_RMSE"] = ["Average_RMSE_INITIAL = " + diff_measures_initial_rmse_res, "Average_RMSE_OPTIMAL = " + diff_measures_optimal_rmse_res ]
1495 __resultats["ResultsSummary_RelativeDifference"] = ["Average_RelativeDifference_INITIAL = " + diff_measures_initial_reldiff_res + "%", "Average_RelativeDifference_OPTIMAL = " + diff_measures_optimal_reldiff_res + "%" ]
1497 measures_input = __resultats ["Measures"]
1498 optimal_simu = __resultats ["OptimalSimulation"]
1500 diff_measures_optimal_rmse = numpy.array([(measures_input[i]-optimal_simu[i])**2 for i in range(len(measures_input))])
1501 diff_measures_optimal_reldiff = numpy.array([abs((measures_input[i]-optimal_simu[i])/measures_input[i]) for i in range(len(measures_input)) if measures_input[i] !=0])
1504 diff_measures_optimal_rmse_res = str(numpy.format_float_scientific(numpy.sqrt(diff_measures_optimal_rmse.mean()), precision = 3))
1505 diff_measures_optimal_reldiff_res = str(numpy.round(diff_measures_optimal_reldiff.mean()*100,decimals =2))
1507 __resultats["ResultsSummary_RMSE"] = ["Average_RMSE_OPTIMAL = " + diff_measures_optimal_rmse_res ]
1508 __resultats["ResultsSummary_RelativeDifference"] = ["Average_RelativeDifference_OPTIMAL = " + diff_measures_optimal_reldiff_res + "%" ]
1513 self.__result["SourceName"],
1514 self.__result["Level"],
1515 self.__result["Format"],
1520 if CreateOptimaldsin:
1521 if self.__model["Format"].upper() in ["DYMOSIM"]:
1522 for etat in range(len(LNames)):
1523 dict_optimal_params = dict(zip(VariablesToCalibrate, __adaocase.get("Analysis")[-1])) #dans la boucle, sinon il se vide à chaque fois que l'on fait itère
1524 if not(ComplexModel):
1525 dir_for_dsin_opti = os.path.abspath(os.path.join(self._get_Name_ModelTmpDir_simple(ModelName),os.pardir))
1526 simu_dir_opti = os.path.abspath(self._get_Name_ModelTmpDir_simple(ModelName))
1528 writer_opti = write_in_dsin.Write_in_dsin(dict_inputs =dict_optimal_params,filedir=simu_dir_opti,dsin_name=str("dsin_" +LNames[etat]+ ".txt"),old_file_name='dsin.txt',new_file_name= str("dsin_" +LNames[etat]+ "_optimal.txt"))
1529 writer_opti.write_in_dsin
1531 writer_opti.write_in_dsin()
1534 os.remove(os.path.join(dir_for_dsin_opti,str("dsin_" +LNames[etat]+ "_optimal.txt")))
1539 os.path.join(simu_dir_opti, str("dsin_" +LNames[etat]+ "_optimal.txt")),
1540 os.path.join(dir_for_dsin_opti, str("dsin_" +LNames[etat]+ "_optimal.txt"))
1543 auto_simul_opti = Automatic_Simulation(
1544 simu_dir=simu_dir_opti,
1545 dsin_name= str("dsin_" +LNames[etat]+ "_optimal.txt"),
1546 dymosim_path = os.path.join(simu_dir_opti, os.pardir) ,
1547 dsres_path = str("dsres_" +LNames[etat]+ "_opti.mat"),
1550 timeout=TimeoutModelExecution,
1551 without_modelicares=True,
1554 auto_simul_opti.single_simulation()
1556 if not(auto_simul_opti.success_code == 0):
1557 print("The dsin.txt with the optimal values does not converge directly, please update the value of iteration variables using Change_working_Point script (consider both calibrated paramaters and boundary conditions in the procedure)")
1559 try: #to handle situations in which there is only one CL (boundary conditions) file
1560 var_int = numpy.ravel(LColumns[:,etat])
1562 var_int = numpy.ravel(LColumns)
1563 dict_inputs_boundary_conditions_optimal_params = dict(zip(LVariablesToChange, var_int))
1564 dict_inputs_boundary_conditions_optimal_params.update(dict_optimal_params)
1567 os.remove(os.path.join(dir_for_dsin_opti,str("dsin_" +LNames[etat]+ "_optimal.txt")))
1570 dir_for_dsin_opti = os.path.abspath(os.path.join(self._get_Name_ModelTmpDir_simple(ModelName),os.pardir))
1572 os.path.join(dir_for_dsin_opti, str("dsin.txt")),
1573 os.path.join(dir_for_dsin_opti, str("dsin_" +LNames[etat]+ "_optimal.txt"))
1575 writer_opti = write_in_dsin.Write_in_dsin(dict_inputs =dict_inputs_boundary_conditions_optimal_params,filedir=dir_for_dsin_opti,dsin_name=str("dsin_" +LNames[etat]+ "_optimal.txt"),old_file_name=str("dsin_" +LNames[etat]+ "_old.txt"),new_file_name= str("dsin_" +LNames[etat]+ "_optimal.txt"))
1576 writer_opti.write_in_dsin
1577 writer_opti.write_in_dsin()
1578 os.remove(os.path.join(dir_for_dsin_opti,str("dsin_" +LNames[etat]+ "_old.txt")))
1579 print("The dsin.txt with the optimal values has not been tested, updating the values of iteration variables might be necessary: this can be done by using Change_working_Point script (consider both calibrated paramaters and boundary conditions in the procedure)")
1581 else: #FMI, OM or GUESS format
1582 print("CreateOptimaldsin option is only valid for DYMOSIM format (it has no effect for the current calculation)")
1584 if self.__verbose: print("")
1587 if not(KeepCalculationFolders):
1588 shutil.rmtree(self._get_Name_ModelTmpDir_simple(ModelName), ignore_errors=True)
1590 if CreateCSVoutputResults:
1591 print("Creation of CSV output results")
1594 optimal_results = __resultats["OptimalSimulation"]
1595 for obs_file_index in range(len(self.__measures["SourceName"])):
1596 obs_filename = self.__measures["SourceName"][obs_file_index] #for the name of the outputresult
1598 time_measures_number = len(OHours[obs_file_index]) #get the number of lines in the obs file
1599 total_measures_number_per_cl = time_measures_number*len(__resultats["NamesOfMeasures"]) #size of the list to be removed from the whole results file
1600 optimal_results_cl = optimal_results[:total_measures_number_per_cl] # output results for the current cl
1601 optimal_results = optimal_results[total_measures_number_per_cl:] # update of the whole results file (we remove te results already treated)
1603 optimal_results_cl_to_csv = optimal_results_cl.reshape((len(__resultats["NamesOfMeasures"]),-1)).transpose() #reshape to get one line per timedate
1604 df_optimal_results_cl_to_csv = pandas.DataFrame(optimal_results_cl_to_csv, columns = __resultats["NamesOfMeasures"]) #creation of the df with the names of measures
1605 df_optimal_results_cl_to_csv.insert(0,"Hour", OHours[obs_file_index]) #add Hour column
1606 df_optimal_results_cl_to_csv.insert(0,"Date", ODates[obs_file_index]) # add Date column
1608 optimal_file_results_name_cl_tmp = "Optimal_simulation_tmp.csv"
1609 df_optimal_results_cl_to_csv.to_csv(optimal_file_results_name_cl_tmp, sep = ";", index = False)
1611 optimal_file_results_name_cl = "Optimal_simulation_"+obs_filename
1612 with open(optimal_file_results_name_cl_tmp, 'r') as infile:
1613 readie=csv.reader(infile, delimiter=';')
1614 with open(optimal_file_results_name_cl, 'wt', newline='') as output:
1615 outwriter=csv.writer(output, delimiter=';')
1616 outwriter.writerow(["#"])
1618 outwriter.writerow(row)
1620 os.remove(optimal_file_results_name_cl_tmp)
1624 if InitialSimulation: #If initial simulation is required as well, the Initial results files are given as well
1625 initial_results = __resultats["InitialSimulation"]
1626 for obs_file_index in range(len(self.__measures["SourceName"])):
1627 obs_filename = self.__measures["SourceName"][obs_file_index] #for the name of the outputresult
1629 time_measures_number = len(OHours[obs_file_index]) #get the number of lines in the obs file
1630 total_measures_number_per_cl = time_measures_number*len(__resultats["NamesOfMeasures"]) #size of the list to be removed from the whole results file
1631 initial_results_cl = initial_results[:total_measures_number_per_cl] # output results for the current cl
1632 initial_results = initial_results[total_measures_number_per_cl:] # update of the whole results file (we remove te results already treated)
1634 initial_results_cl_to_csv = initial_results_cl.reshape((len(__resultats["NamesOfMeasures"]),-1)).transpose() #reshape to get one line per timedate
1635 df_initial_results_cl_to_csv = pandas.DataFrame(initial_results_cl_to_csv, columns = __resultats["NamesOfMeasures"]) #creation of the df with the names of measures
1636 df_initial_results_cl_to_csv.insert(0,"Hour", OHours[obs_file_index]) #add Hour column
1637 df_initial_results_cl_to_csv.insert(0,"Date", ODates[obs_file_index]) # add Date column
1639 initial_file_results_name_cl_tmp = "Initial_simulation_tmp.csv"
1640 df_initial_results_cl_to_csv.to_csv(initial_file_results_name_cl_tmp, sep = ";", index = False)
1642 initial_file_results_name_cl = "Initial_simulation_"+obs_filename
1643 with open(initial_file_results_name_cl_tmp, 'r') as infile:
1644 readie=csv.reader(infile, delimiter=';')
1645 with open(initial_file_results_name_cl, 'wt', newline='') as output:
1646 outwriter=csv.writer(output, delimiter=';')
1647 outwriter.writerow(["#"])
1649 outwriter.writerow(row)
1651 os.remove(initial_file_results_name_cl_tmp)
1658 def verify_gradient(self):
1659 raise NotImplementedError("Not yet implemented")
1661 # ==============================================================================
1662 def _readLink(__filename=None, __colnames=None, __indexname="Variable", __format="Guess"):
1664 Read file of link between model and measures
1668 - Names of the columns to read, known by their variable header
1669 - Name of the column containing the variables to adapt
1670 - Format of the data (CSV, TSV... or can be guessed)
1673 if __colnames is None:
1675 if __indexname is None: __indexname="Variable"
1677 try: #Solution provisoire pou gérer les cas sans CL précisées par l'utilisateur
1678 with ImportFromFile(__filename, __colnames, __indexname, __format, False) as reading:
1679 colnames, columns, indexname, variablestochange = reading.getvalue()
1681 colnames =__colnames
1683 variablestochange =[]
1684 print("/!\\ Boundary conditions indicated in ", __filename, " are not considered for calculation. Indicate at least two boundary conditions so that they are considered. /!\\ ")
1686 variablestochange = tuple([v.strip() for v in variablestochange])
1687 return (colnames, columns, variablestochange)
1689 # ==============================================================================
1690 def _readMultiMeasures(__filenames=None, __colnames=None, __format="Guess"):
1692 Read files of measures, using only some named variables by columns
1696 - Names of the columns to read, known by their variable header
1697 - Format of the data (CSV, TSV... or can be guessed)
1701 for __filename in __filenames:
1702 colnames, columns = _readMeasures(__filename, __colnames, __format)
1703 MultiMeasures.append( columns )
1705 return (colnames, MultiMeasures)
1707 # ==============================================================================
1708 def _readMeasures(__filename=None, __colnames=None, __format="Guess"):
1710 Read files of measures, using only some named variables by columns
1714 - Names of the columns to read, known by their variable header
1715 - Format of the data (CSV, TSV... or can be guessed)
1718 with ImportFromFile(__filename, __colnames, None, __format) as reading:
1719 colnames, columns, indexname, index = reading.getvalue()
1721 return (colnames, columns)
1723 # ==============================================================================
1724 def _readBackground(__filename="dsin.txt", __varnames=None, __backgroundformat="Guess"):
1726 Read background in model definition. The priority is "USER", then "ADAO",
1727 then "DSIN", and the default one correspond to "DSIN".
1731 - Names of the variables to be found in the model data
1732 - Format of the data (can be guessed)
1735 if __backgroundformat.upper() in ["GUESS", "DSIN"]:
1737 if __filename is not None and os.path.isfile(__filename):
1738 __model_dir = os.path.dirname(__filename)
1739 __model_fil = os.path.basename(__filename)
1740 elif __filename is not None and os.path.isdir(__filename) and os.path.isfile(os.path.join(__filename, "dsin.txt")):
1741 __model_dir = os.path.abspath(__filename)
1742 __model_fil = "dsin.txt"
1744 raise ValueError('No such file or directory: %s'%str(__filename))
1745 __bgfile = os.path.abspath(os.path.join(__model_dir, __model_fil))
1746 if (len(__bgfile)<7) and (__bgfile[-8:].lower() != "dsin.txt" ):
1747 raise ValueError("Model definition file \"dsin.txt\" can not be found.")
1748 if __varnames is None:
1749 raise ValueError("Model variables to be read has to be given")
1750 if __backgroundformat.upper() == "ADAO":
1752 if __filename is None or not os.path.isfile(__filename):
1753 raise ValueError('No such file or directory: %s'%str(__filename))
1754 __bgfile = os.path.abspath(__filename)
1755 if not( len(__bgfile)>4 ):
1756 raise ValueError("User file name \"%s\" is too short and seems not to be valid."%__bgfile)
1757 if __backgroundformat.upper() == "USER":
1759 if __filename is None or not os.path.isfile(__filename):
1760 raise ValueError('No such file or directory: %s'%str(__filename))
1761 __bgfile = os.path.abspath(__filename)
1762 if not( len(__bgfile)>5 ):
1763 raise ValueError("User file name \"%s\" is too short and seems not to be valid."%__bgfile)
1764 if __format not in ["DSIN", "ADAO", "USER"]:
1765 raise ValueError("Background initial values asked format \"%s\" is not valid."%__format)
1767 #---------------------------------------------
1768 if __format == "DSIN":
1769 if __varnames is None:
1770 raise ValueError("Names of variables to read has to be given")
1772 __names, __background, __bounds = [], [], []
1773 with open(__bgfile, 'r') as fid:
1774 __content = fid.readlines()
1775 for v in tuple(__varnames):
1776 if _verify_existence_of_variable_in_dsin(v, __content):
1777 __value = _read_existing_variable_in_dsin(v, __content)
1779 __background.append(__value[0])
1780 __bounds.append(__value[1:3])
1781 # print "%s = "%v, __value
1782 __names = tuple(__names)
1783 __background = numpy.ravel(__background)
1784 __bounds = tuple(__bounds)
1785 #---------------------------------------------
1786 elif __format.upper() == "ADAO":
1787 with open(__bgfile, 'r') as fid:
1790 __background = locals().get("Background", None)
1791 if __background is None:
1792 raise ValueError("Background is not defined")
1794 __background = numpy.ravel(__background)
1795 __Params = locals().get("Parameters", {})
1796 if "Bounds" not in __Params:
1797 raise ValueError("Bounds can not be find in file \"%s\""%str(__filename))
1799 __bounds = tuple(__Params["Bounds"])
1800 if len(__bounds) != len(__background):
1801 raise ValueError("Initial background length does not match bounds number")
1803 if __varnames is None:
1806 __names = tuple(__varnames)
1807 if len(__names) != len(__background):
1808 raise ValueError("Initial background length does not match names number")
1809 #---------------------------------------------
1810 elif __format.upper() == "USER":
1811 __names, __background, __bounds = ImportScalarLinesFromFile(__bgfile).getvalue(__varnames)
1812 #---------------------------------------------
1813 for i,b in enumerate(__bounds) :
1814 if b[0] is not None and b[1] is not None and not (b[0] < b[1]) :
1815 raise ValueError(f'Inconsistent boundaries values for "{__names[i]}": {b[0]} not < {b[1]}')
1817 return (__names, __background, __bounds)
1819 def _verify_existence_of_variable_in_dsin(__variable, __content, __number = 2):
1820 "Verify if the variable exist in the model file content"
1821 if "".join(__content).count(__variable) >= __number:
1826 def _read_existing_variable_in_dsin(__variable, __content):
1827 "Return the value of the real variable"
1828 __previousline, __currentline = "", ""
1830 initialValuePart = False
1831 for line in __content:
1832 if not initialValuePart and "initialValue" not in line: continue
1833 if not initialValuePart and "initialValue" in line: initialValuePart = True
1834 __previousline = __currentline
1835 __currentline = line
1836 if __variable in __currentline and "#" in __currentline and "#" not in __previousline:
1837 # Informations ecrites sur deux lignes successives
1838 vini, value, vmin, vmax = __previousline.split()
1839 # vcat, vtype, diez, vnam = __currentline.split()
1840 value = float(value)
1841 if float(vmin) >= float(vmax):
1842 vmin, vmax = None, None
1844 vmin, vmax = float(vmin), float(vmax)
1845 return (value, vmin, vmax)
1846 elif __variable in __currentline and "#" in __currentline and "#" in __previousline:
1847 # Informations ecrites sur une seule ligne
1848 vini, value, vmin, vmax, reste = __previousline.split(maxsplit=5)
1849 value = float(value)
1850 if float(vmin) >= float(vmax):
1851 vmin, vmax = None, None
1853 vmin, vmax = float(vmin), float(vmax)
1854 return (value, vmin, vmax)
1856 raise ValueError("Value of variable %s not found in the file content"%__variable)
1858 # ==============================================================================
1859 def _readMethod(__filename="parameters.py", __format="ADAO"):
1861 Read data assimilation method parameters
1865 - Format of the data (can be guessed)
1867 if not os.path.isfile(__filename):
1868 raise ValueError('No such file or directory: %s'%str(__filename))
1870 #---------------------------------------------
1871 if __format.upper() in ["GUESS", "ADAO"]:
1872 with open(__filename, 'r') as fid:
1874 __Algo = locals().get("Algorithm", "3DVAR")
1875 __Params = locals().get("Parameters", {})
1876 __CovB = locals().get("BackgroundError", 1.)
1877 __CovR = locals().get("ObservationError", 1.)
1878 # __Xb = locals().get("Background", __background)
1879 # if not isinstance(__Params, dict): __Params = {"Bounds":__bounds}
1880 # if "Bounds" not in __Params: __Params["Bounds"] = __bounds
1881 #---------------------------------------------
1883 __Algo, __Params, __CovB, __CovR = "3DVAR", {}, 1., 1.
1884 #---------------------------------------------
1886 # raise ValueError("Background is not defined")
1888 # __Xb = numpy.ravel(__Xb)
1890 return (__Algo, __Params, __CovB, __CovR)
1892 # ==============================================================================
1893 def _saveResults(__resultats, __filename=None, __level=None, __format="Guess", __verbose=False):
1895 Save results relatively to the level required
1899 - Level of saving : final output only or intermediary with final output
1900 - Format of the data
1902 if __filename is None:
1903 raise ValueError('A file to save results has to be named, please specify one')
1904 if __format.upper() == "GUESS":
1905 if __filename.split(".")[-1].lower() == "txt":
1907 elif __filename.split(".")[-1].lower() == "py":
1910 raise ValueError("Can not guess the file format of \"%s\", please specify the good one"%__filename)
1912 __format = str(__format).upper()
1913 if __format not in ["TXT", "PY"]:
1914 raise ValueError("Result file format \"%s\" is not valid."%__format)
1915 if __format == "PY" and os.path.splitext(__filename)[1] != ".py":
1916 __filename = os.path.splitext(__filename)[0] + ".py"
1918 #---------------------------------------------
1919 if __format.upper() == "TXT":
1920 output = ["# Final results",]
1921 keys = list(__resultats.keys())
1924 output.append("%18s = %s"%(k,tuple(__resultats[k])))
1926 with open(__filename, 'w') as fid:
1927 fid.write( "\n".join(output) )
1929 #---------------------------------------------
1930 elif __format.upper() == "PY":
1931 output = ["# Final results",]
1932 keys = list(__resultats.keys())
1935 output.append("%s = %s"%(k,tuple(__resultats[k])))
1937 with open(__filename, 'w') as fid:
1938 fid.write( "\n".join(output) )
1940 #---------------------------------------------
1941 elif __format.upper() == "CSV":
1942 raise NotImplementedError("Not yet implemented")
1943 #---------------------------------------------
1944 elif __format.upper() == "DICT":
1945 raise NotImplementedError("Not yet implemented")
1946 #---------------------------------------------
1947 if __verbose: # Format TXT
1948 output = ["# Final results",]
1949 keys = list(__resultats.keys())
1952 output.append("%18s = %s"%(k,tuple(__resultats[k])))
1954 print( "\n".join(output) )
1955 #---------------------------------------------
1958 def _secure_copy_textfile(__original_file, __destination_file):
1959 "Copy the file guranteeing that it is copied"
1966 shutil.copystat( # commande pour assurer que le fichier soit bien copié
1970 os.path.getsize(__destination_file)# commande bis pour assurer que le fichier soit bien copié
1972 def file_len(fname):
1974 with open(fname) as f:
1979 while file_len(__original_file)!= file_len(__destination_file):
1985 # ==============================================================================
1986 class Python_Simulation(object):
1987 def __init__(self, simu_dir=".", val_inputs=(), logfile=True, timeout=60, verbose=True):
1988 self.__simu_dir = os.path.abspath(simu_dir)
1989 self.__inputs = numpy.ravel(val_inputs)
1992 print(" [VERBOSE] Input values %s"%self.__inputs)
1994 def single_simulation(self, __inputs = None):
1995 with open(os.path.join(self.__simu_dir,"pythonsim.exe"), 'r') as fid:
1997 __directoperator = locals().get("DirectOperator", None)
1998 if __inputs is None: __inputs = self.__inputs
1999 return __directoperator(__inputs)
2001 # ==============================================================================
2002 class VersionsLogicielles (object):
2003 "Modules version information"
2005 print(" [VERBOSE] System configuration:")
2006 print(" - Python...:",sys.version.split()[0])
2007 print(" - Numpy....:",numpy.version.version)
2008 print(" - Scipy....:",scipy.version.version)
2009 print(" - ADAO.....:",adao.version)
2012 def TOP_LEVEL_exefmuMultiobs_ref( x_values_matrix , VariablesToCalibrate=None, OutputVariables=None, ref_simudir = None, ModelName = None ):
2013 "Appel du modèle en format FMU et restitution des résultats"
2014 x_values = list(numpy.ravel(x_values_matrix))
2015 x_inputs=dict(zip(VariablesToCalibrate,x_values))
2017 fmuname = os.path.basename(ModelName) #in case ModelName is given as a path (works also if it is a file name)
2018 fmu = os.path.join(ref_simudir,fmuname)
2020 reader = simulate_fmu(fmu, output = OutputVariables, start_values = x_inputs)
2021 y_values = [reader[y_name][-1] for y_name in OutputVariables]
2022 y_values_matrix = numpy.asarray(y_values)
2024 return y_values_matrix
2026 def TOP_LEVEL_exefmuMultiobs( x_values_matrix , VariablesToCalibrate=None, OutputVariables=None, LNames = None, LColumns = None, LVariablesToChange=None, ref_simudir = None, ModelName = None, List_Multideltatime = None, AdvancedDebugModel = None, TimeoutModelExecution = None, FMUInput = None):
2027 "Appel du modèle en format FMU et restitution des résultats"
2029 if VariablesToCalibrate is None:
2030 raise ValueError("VariablesToCalibrate")
2032 if OutputVariables is None:
2033 raise ValueError("OutputVariables")
2036 raise ValueError("LNames")
2038 if LColumns is None:
2039 raise ValueError("LColumns")
2041 if LVariablesToChange is None:
2042 raise ValueError("LVariablesToChange")
2044 if ref_simudir is None:
2045 raise ValueError("ref_simudir")
2047 if List_Multideltatime is None:
2048 raise ValueError("Problem defining simulation output results")
2050 if AdvancedDebugModel is None:
2051 raise ValueError("AdvancedDebugModel")
2053 if TimeoutModelExecution is None:
2054 raise ValueError("TimeoutModelExecution")
2056 x_values = list(numpy.ravel(x_values_matrix))
2057 x_inputs=dict(zip(VariablesToCalibrate,x_values))
2059 if os.path.isfile(ModelName):
2060 fmuname = os.path.basename(ModelName)
2061 elif os.path.isdir(ModelName):
2062 fmuname = [files for files in os.listdir(os.path.abspath(ModelName)) if files[-4:] =='.fmu'][0]
2064 fmu = os.path.join(ref_simudir,fmuname)
2066 multi_y_values_matrix = []
2068 for etat in range(len(LNames)):
2071 try: #to handle situations in which there is only one CL (boundary conditions) file
2072 var_int = numpy.ravel(LColumns[:,etat])
2074 var_int = numpy.ravel(LColumns)
2076 dict_inputs = dict(zip(LVariablesToChange, var_int))
2078 dict_inputs.update(x_inputs)
2081 fmu_inputs = FMUInput[LNames[etat]]
2082 timestep = fmu_inputs[1][0] - fmu_inputs[0][0] #Assuming constant timestep
2083 if AdvancedDebugModel: print(f'The timestep for {LNames[etat]} is {timestep} seconds')
2089 stoptime_fmu = List_Multideltatime[etat][-1]
2095 if AdvancedDebugModel == True:
2096 old_stdout = sys.stdout
2097 new_stdout = io.StringIO()
2098 sys.stdout = new_stdout
2100 start_time_simulation_fmi = time.time()#timeout manangement since fmpy does not raise an error for this, it just ends the simulations and continues
2103 reader = simulate_fmu(fmu, output = OutputVariables, start_values = dict_inputs, debug_logging = True, timeout = TimeoutModelExecution, input = fmu_inputs, stop_time= stoptime_fmu, output_interval= timestep)
2104 except FMICallException as e:
2105 output = new_stdout.getvalue()
2106 sys.stdout = old_stdout
2107 print('[ERROR] Failed simulation with the following input:\n',dict_inputs)
2110 elapsed_time_simulation_fmi = time.time() - start_time_simulation_fmi
2111 if elapsed_time_simulation_fmi > TimeoutModelExecution:
2112 raise TimeoutError("Timeout for simulation reached, please increase it in order to be able to simulate your model and/or check if your model is correct (use TimeoutModelExecution option in configuration.py file)")
2114 output = new_stdout.getvalue()
2115 sys.stdout = old_stdout
2118 log_file=os.path.join(dir_run,
2119 'log_debug_model.txt')
2120 try: #try toremove the previous file
2125 f=open(log_file,'a')
2130 start_time_simulation_fmi = time.time()
2133 reader = simulate_fmu(fmu, output = OutputVariables, start_values = dict_inputs, timeout = TimeoutModelExecution, input = fmu_inputs, stop_time= stoptime_fmu, output_interval= timestep)
2134 except FMICallException as e:
2135 print('[ERROR] Failed simulation with the following input:\n',dict_inputs)
2138 elapsed_time_simulation_fmi = time.time() - start_time_simulation_fmi
2139 if elapsed_time_simulation_fmi > TimeoutModelExecution:
2140 raise TimeoutError("Timeout for simulation reached, please increase it in order to be able to simulate your model and/or check if your model is correct (use TimeoutModelExecution option in configuration.py file)")
2142 y_whole = [reader[y_name] for y_name in OutputVariables]
2144 for y_ind in range(len(OutputVariables)):
2145 y_ind_whole_values = y_whole[y_ind]
2146 y_ind_whole_time = reader['time'] #standard output of fmu
2148 if len(List_Multideltatime[etat])>1:
2149 index_y_values = [find_nearest_index(y_ind_whole_time,time)[0] for time in List_Multideltatime[etat]]
2151 index_y_values = [-1] #we only take the last point if one measure point
2152 y_ind_values = [y_ind_whole_values[i] for i in index_y_values]
2153 y_values = y_values + y_ind_values
2155 y_values_matrix = y_values #pbs in results management
2156 multi_y_values_matrix = multi_y_values_matrix + y_values_matrix
2158 y_values_matrix_def = numpy.ravel(numpy.array(multi_y_values_matrix))
2159 return y_values_matrix_def
2161 def TOP_LEVEL_exeOpenModelicaMultiobs( x_values_matrix , KeepCalculationFolders = None, VariablesToCalibrate=None, OutputVariables=None, LNames = None, LColumns = None, LVariablesToChange=None, ref_simudir = None, ModelName = None, List_Multideltatime = None, Linux = None, AdvancedDebugModel = None, TimeoutModelExecution = None):
2162 "Appel du modèle en format OpenModelica et restitution des résultats"
2164 if VariablesToCalibrate is None:
2165 raise ValueError("VariablesToCalibrate")
2167 if OutputVariables is None:
2168 raise ValueError("OutputVariables")
2171 raise ValueError("LNames")
2173 if LColumns is None:
2174 raise ValueError("LColumns")
2176 if LVariablesToChange is None:
2177 raise ValueError("LVariablesToChange")
2179 if ref_simudir is None:
2180 raise ValueError("ref_simudir")
2182 if KeepCalculationFolders is None:
2183 raise ValueError("KeepCalculationFolders")
2186 raise ValueError("Linux")
2188 if List_Multideltatime is None:
2189 raise ValueError("Problem defining simulation output results")
2191 if AdvancedDebugModel is None:
2192 raise ValueError("AdvancedDebugModel")
2194 if TimeoutModelExecution is None:
2195 raise ValueError("TimeoutModelExecution")
2197 x_values = list(numpy.ravel(x_values_matrix))
2198 x_inputs=dict(zip(VariablesToCalibrate,x_values))
2200 if os.path.isfile(ModelName):
2201 om_xml = os.path.basename(ModelName)
2202 elif os.path.isdir(ModelName):
2203 om_xml = [files for files in os.listdir(os.path.abspath(ModelName)) if files[-4:] =='.xml'][0]
2205 om_xml = os.path.abspath(os.path.join(ref_simudir,om_xml))
2207 multi_y_values_matrix = []
2209 for etat in range(len(LNames)):
2212 try: #to handle situations in which there is only one CL (boundary conditions) file
2213 var_int = numpy.ravel(LColumns[:,etat])
2215 var_int = numpy.ravel(LColumns)
2217 dict_inputs = dict(zip(LVariablesToChange, var_int))
2220 dict_inputs.update(x_inputs)
2223 results_dir = tempfile.mkdtemp( dir=os.path.dirname(om_xml) ) #to handle parallel computation
2225 results_file_name = os.path.join(results_dir,os.path.basename(om_xml)[:-9] + '_' + LNames[etat].replace('.','_') + ".mat")
2227 OM_execution = run_OM_model(om_xml, dict_inputs = dict_inputs, result_file_name = results_file_name , Linux = Linux, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution)
2230 reader = Reader(results_file_name,'dymola') #dymola even if it is OpenModelica
2232 raise ValueError("Simulation cannot be performed: reduce the number of parameters to calibrate and/or the range in which their optimal value should be found (or modify and simplify the model to make it easier to simulate)" )
2234 y_whole = [reader.values(y_name) for y_name in OutputVariables]
2236 for y_ind in range(len(OutputVariables)):
2237 y_ind_whole_values = y_whole[y_ind][1]
2238 y_ind_whole_time = y_whole[y_ind][0]
2240 if len(List_Multideltatime[etat])>1:
2241 index_y_values = [find_nearest_index(y_ind_whole_time,time)[0] for time in List_Multideltatime[etat]]
2243 index_y_values = [-1] #we only take the last point if one measure point
2245 y_ind_values = [y_ind_whole_values[i] for i in index_y_values]
2246 y_values = y_values + y_ind_values
2247 # y_values_matrix = numpy.asarray(y_values)
2248 y_values_matrix = y_values #pbs in results management
2249 multi_y_values_matrix = multi_y_values_matrix + y_values_matrix
2251 if not KeepCalculationFolders:
2252 shutil.rmtree(results_dir, ignore_errors=True)
2254 y_values_matrix_def = numpy.ravel(numpy.array(multi_y_values_matrix))
2256 return y_values_matrix_def
2258 def TOPLEVEL_exedymosimMultiobs_simple( x_values_matrix , VariablesToCalibrate=None, OutputVariables=None, LNames = None, ref_simudir = None, KeepCalculationFolders = None, Linux = None, List_Multideltatime = None, AdvancedDebugModel = None, TimeoutModelExecution = None):
2259 "Appel du modèle en format DYMOSIM et restitution des résultats"
2262 if VariablesToCalibrate is None:
2263 raise ValueError("VariablesToCalibrate")
2265 if OutputVariables is None:
2266 raise ValueError("OutputVariables")
2269 raise ValueError("LNames")
2271 if ref_simudir is None:
2272 raise ValueError("ref_simudir")
2274 if KeepCalculationFolders is None:
2275 raise ValueError("KeepCalculationFolders")
2278 raise ValueError("Linux")
2280 if List_Multideltatime is None:
2281 raise ValueError("Problem defining simulation output results")
2283 if AdvancedDebugModel is None:
2284 raise ValueError("AdvancedDebugModel")
2286 if TimeoutModelExecution is None:
2287 raise ValueError("TimeoutModelExecution")
2289 # x_values = list(numpy.ravel(x_values_matrix) * Background)
2290 x_values = list(numpy.ravel(x_values_matrix))
2291 x_inputs=dict(zip(VariablesToCalibrate,x_values))
2292 # if Verbose: print(" Simulation for %s"%numpy.ravel(x_values))
2295 multi_y_values_matrix = []
2297 for etat in range(len(LNames)):
2300 simudir = tempfile.mkdtemp( dir=ref_simudir )
2303 dict_inputs.update( x_inputs )
2306 os.path.join(ref_simudir, str("dsin_" +LNames[etat]+ ".txt")),
2307 os.path.join(simudir, str("dsin_" +LNames[etat]+ ".txt"))
2309 _secure_copy_textfile(
2310 os.path.join(ref_simudir, str("dsin_" +LNames[etat]+ ".txt")),
2311 os.path.join(simudir, str("dsin_" +LNames[etat]+ ".txt"))
2314 auto_simul = Automatic_Simulation(
2316 dymosim_path = os.path.join(ref_simudir,os.pardir) ,
2317 dsin_name = str("dsin_" +LNames[etat]+ ".txt"),
2318 dsres_path = os.path.join(simudir, str("dsres" +LNames[etat]+ ".mat")),
2319 dict_inputs = dict_inputs,
2321 timeout = TimeoutModelExecution,
2322 without_modelicares = True,
2325 auto_simul.single_simulation()
2327 if AdvancedDebugModel:
2328 dslog_file = os.path.join(simudir, 'dslog.txt')
2329 debug_model_file = os.path.join(ref_simudir,'log_debug_model.txt')
2331 shutil.copy(dslog_file,debug_model_file)
2335 if auto_simul.success_code == 2 :
2336 raise ValueError("The simulation falis after initialization, please check the log file and your model in order to be sure that the whole simulation can be performed (specially for dynamic models)")
2338 if auto_simul.success_code == 3 :
2339 raise ValueError("The simulation did not reach the final time, probably due to a lack of time, please increase the time for model execution by modifying the TimeoutModelExecution value in configuration.py")
2341 if auto_simul.success_code == 0:
2342 pass #means that everything OK, we keep going
2344 raise ValueError("The NeedInitVariables iteration should be set to True, convergence issues found for %s" % str(LNames[etat])) #means this value is 1--> not correct initialization
2346 reader = Reader(os.path.join(simudir, str("dsres" +LNames[etat]+ ".mat")),'dymola')
2348 y_whole = [reader.values(y_name) for y_name in OutputVariables]
2350 for y_ind in range(len(OutputVariables)):
2351 y_ind_whole_values = y_whole[y_ind][1]
2352 y_ind_whole_time = y_whole[y_ind][0]
2354 if len(List_Multideltatime[etat])>1:
2355 index_y_values = [find_nearest_index(y_ind_whole_time,time)[0] for time in List_Multideltatime[etat]]
2357 index_y_values = [-1] #we only take the last point if one measure point
2359 y_ind_values = [y_ind_whole_values[i] for i in index_y_values]
2360 y_values = y_values + y_ind_values
2361 # y_values_matrix = numpy.asarray(y_values)
2362 y_values_matrix = y_values #pbs in results management
2363 multi_y_values_matrix = multi_y_values_matrix + y_values_matrix
2365 if not KeepCalculationFolders:
2366 shutil.rmtree(simudir, ignore_errors=True)
2368 y_values_matrix_def = numpy.ravel(numpy.array(multi_y_values_matrix))
2369 return y_values_matrix_def
2371 def TOP_LEVEL_exedymosimMultiobs( x_values_matrix , VariablesToCalibrate=None, OutputVariables=None, LNames = None, ModelFormat = None, KeepCalculationFolders = None, Verbose = None, dict_dsin_paths = None, Linux = None, List_Multideltatime = None, AdvancedDebugModel = None, TimeoutModelExecution = None): #2ème argument à ne pas modifier
2372 "Appel du modèle en format DYMOSIM et restitution des résultats"
2374 if VariablesToCalibrate is None:
2375 raise ValueError("VariablesToCalibrate")
2377 if OutputVariables is None:
2378 raise ValueError("OutputVariables")
2381 raise ValueError("LNames")
2383 if ModelFormat is None:
2384 raise ValueError("ModelFormat")
2386 if KeepCalculationFolders is None:
2387 raise ValueError("KeepCalculationFolders")
2390 raise ValueError("Verbose")
2392 if dict_dsin_paths is None:
2393 raise ValueError("dict_dsin_paths")
2396 raise ValueError("Linux")
2398 if List_Multideltatime is None:
2399 raise ValueError("Problem defining simulation output results")
2401 if AdvancedDebugModel is None:
2402 raise ValueError("AdvancedDebugModel")
2404 if TimeoutModelExecution is None:
2405 raise ValueError("TimeoutModelExecution")
2407 # x_values = list(numpy.ravel(x_values_matrix) * Background)
2408 x_values = list(numpy.ravel(x_values_matrix))
2409 x_inputs = dict(zip(VariablesToCalibrate,x_values))
2410 x_inputs_for_CWP = {}
2411 # if Verbose: print(" Simulation for %s"%numpy.ravel(x_values))
2413 multi_y_values_matrix = []
2417 for etat in range(len(LNames)):
2420 ref_simudir = Calibration._setModelTmpDir( None,dict_dsin_paths[LNames[etat]], ModelFormat, 'dsin.txt', os.path.join(os.path.pardir,os.path.pardir), LNames[etat])
2421 simudir = ref_simudir
2423 dict_inputs.update(x_inputs)
2426 auto_simul = Automatic_Simulation(
2428 dymosim_path = os.path.join(simudir, os.pardir, os.pardir) ,
2429 dsres_path = os.path.join(simudir, str("dsres" +LNames[etat]+ ".mat")),
2430 dict_inputs=dict_inputs,
2432 timeout=TimeoutModelExecution,
2433 without_modelicares=True,
2436 auto_simul.single_simulation()
2438 if AdvancedDebugModel:
2439 dslog_file = os.path.join(simudir, 'dslog.txt')
2440 debug_model_file = os.path.join(simudir, os.pardir, 'log_debug_model.txt')
2442 shutil.copy(dslog_file,debug_model_file)
2446 if auto_simul.success_code == 2 :
2447 raise ValueError("The simulation falis after initialization, please check the log file and your model in order to be sure that the whole simulation can be performed (specially for dynamic models)")
2449 if auto_simul.success_code == 3 :
2450 raise ValueError("The simulation did not reach the final time, probably due to a lack of time, please increase the time for model execution by modifying the TimeoutModelExecution value in configuration.py")
2453 if auto_simul.success_code == 0:
2454 reader = Reader(os.path.join(simudir, str("dsres" +LNames[etat]+ ".mat")),'dymola')
2458 path_around_simu = os.path.join(simudir,os.path.pardir,os.path.pardir)
2459 around_simu = Around_Simulation(dymola_version='2012',
2460 curdir=path_around_simu,
2461 source_list_iter_var = 'from_dymtranslog',
2462 copy_from_dym_trans_log= os.path.join(path_around_simu,'ini.txt'),
2463 mat = os.path.join(os.path.dirname(dict_dsin_paths[LNames[etat]]), str("REF_for_CWP" + ".mat")),
2464 iter_var_values_options={'source':'from_mat','moment':'initial'},
2466 around_simu.set_list_iter_var()
2467 around_simu.set_dict_iter_var()
2469 reader_CWP = Reader(os.path.join(os.path.dirname(dict_dsin_paths[LNames[etat]]), str("REF_for_CWP" + ".mat")),'dymola')
2470 x_values_intemediary = [reader_CWP.values(x_name)[1][-1] for x_name in VariablesToCalibrate]
2471 x_values_from_last_calculation = numpy.asarray(x_values_intemediary)
2472 x_inputs_from_last_calculation = dict(zip(VariablesToCalibrate,x_values_from_last_calculation))
2475 for var_calib in x_inputs_from_last_calculation.keys():
2476 x_inputs_for_CWP[var_calib] = (x_inputs_from_last_calculation[var_calib], x_inputs[var_calib])
2478 dict_var_to_fix2 = Dict_Var_To_Fix(option='automatic',
2479 dict_auto_var_to_fix=x_inputs_for_CWP)
2481 dict_var_to_fix2.set_dict_var_to_fix()
2484 LOG_FILENAME = os.path.join(simudir,'change_working_point.log')
2485 for handler in logging.root.handlers[:]:
2486 logging.root.removeHandler(handler)
2487 logging.basicConfig(filename=LOG_FILENAME,level=logging.DEBUG,filemode='w')
2489 working_point_modif = Working_Point_Modification(main_dir = simudir,
2490 simu_material_dir = simudir,
2491 dymosim_path = os.path.join(simudir, os.pardir, os.pardir),
2493 store_res_dir = 'RES',
2494 dict_var_to_fix = dict_var_to_fix2.dict_var_to_fix,
2495 around_simulation0 = around_simu,
2496 var_to_follow_path= os.path.join(simudir,'var_to_follow.csv'),
2497 gen_scripts_ini= False,
2498 nit_max= 1000000000000000,
2499 min_step_val = 0.00000000000005,
2500 timeout = TimeoutModelExecution,
2503 working_point_modif.create_working_directory()
2504 working_point_modif.working_point_modification(skip_reference_simulation=True)
2506 if AdvancedDebugModel:
2507 dslog_file = os.path.join(simudir,'SIMUS', 'dslog.txt')
2508 debug_model_file = os.path.join(simudir, os.pardir, 'log_debug_model.txt')
2510 shutil.copy(dslog_file,debug_model_file)
2514 reader = Reader(os.path.join(simudir, 'RES','1.0.mat'),'dymola')
2516 raise ValueError("Simulation cannot be performed: reduce the number of parameters to calibrate and/or the range in which their optimal value should be found (or modify and simplify the model to make it easier to simulate)" )
2519 y_whole = [reader.values(y_name) for y_name in OutputVariables]
2521 for y_ind in range(len(OutputVariables)):
2522 y_ind_whole_values = y_whole[y_ind][1]
2523 y_ind_whole_time = y_whole[y_ind][0]
2524 if len(List_Multideltatime[etat])>1:
2525 index_y_values = [find_nearest_index(y_ind_whole_time,time)[0] for time in List_Multideltatime[etat]]
2527 index_y_values = [-1] #we only take the last point if one measure point
2529 y_ind_values = [y_ind_whole_values[i] for i in index_y_values]
2530 y_values = y_values + y_ind_values
2531 # y_values_matrix = numpy.asarray(y_values)
2532 y_values_matrix = y_values #pbs in results management
2533 multi_y_values_matrix = multi_y_values_matrix + y_values_matrix
2535 if not KeepCalculationFolders:
2536 shutil.rmtree(simudir, ignore_errors=True)
2537 y_values_matrix_def = numpy.ravel(numpy.array(multi_y_values_matrix))
2538 return y_values_matrix_def
2540 def run_OM_model(xml_file, #req arg, file path to _init.xml file
2541 result_file_name = None, #file path for storing results
2542 dict_inputs = None, #optionnal argument : for overriding paramters
2544 AdvancedDebugModel = None,
2545 TimeoutModelExecution = None):
2547 Command=[] #initialize command with list of argument
2549 xml_file = os.path.abspath(xml_file)
2551 #set main command to be runned
2553 OM_binary = xml_file.replace('_init.xml','') # not tested
2555 OM_binary = xml_file.replace('_init.xml','.exe')
2556 Command.append(OM_binary)
2558 #Get base path for binaries and input
2561 inputPath='-inputPath='+os.path.dirname(xml_file)
2562 Command.append(inputPath)
2564 #Generate override command
2565 if dict_inputs !=None:
2566 Override='-override='
2567 for elt in dict_inputs.items():
2568 Override=Override+elt[0]+'='+str(elt[1])+','
2569 Override=Override.rstrip(',') #remove last ',' if any
2570 Command.append(Override)
2572 #Generate name of result file
2573 if result_file_name != None:
2574 result='-r='+result_file_name
2575 Command.append(result)
2577 result='-r='+OM_binary+'.mat' #by default result is stored next to _init and bin file
2578 Command.append(result)
2580 result='-w' #by default result is stored next to _init and bin file
2581 Command.append(result)
2583 result='-lv='+'LOG_STATS,LOG_NLS_V' #by default result is stored next to _init and bin file
2584 Command.append(result)
2586 inputPath='-outputPath='+os.path.dirname(xml_file)
2587 Command.append(inputPath)
2590 proc = subprocess.Popen(Command, #Command in the form of a text list
2591 stderr=subprocess.STDOUT,
2592 stdout=subprocess.PIPE,
2593 universal_newlines=True)
2594 # for line in proc.stdout.readlines():
2597 if AdvancedDebugModel == True:
2598 dir_run=os.path.join(os.path.dirname(result_file_name),os.pardir)
2599 log_file=os.path.join(dir_run,
2600 'log_debug_model.txt')
2607 f=open(log_file,'a')
2608 for line in proc.stdout.readlines():
2611 proc.communicate(timeout = TimeoutModelExecution)
2612 except subprocess.TimeoutExpired:
2613 raise ValueError("Timeout for simulation reached, please increase it in order to be able to simulate your model and/or check if your model is correct (use TimeoutModelExecution option in configuration.py file)")
2617 proc.communicate(timeout = TimeoutModelExecution)
2619 raise ValueError("Timeout for simulation reached, please increase it in order to be able to simulate your model and/or check if your model is correct (use TimeoutModelExecution option in configuration.py file)")
2621 # proc.communicate()
2627 # if proc.returncode !=0:
2628 # raise ValueError("Simulation not ended")
2631 # sys.stdout.flush()
2632 # sys.stderr.flush()
2634 # dir_run=os.path.dirname(result_file_name)
2635 # log_file=os.path.join(dir_run,
2638 # proc.wait(timeout=50)
2640 # f=open(log_file,'a')
2641 # for line in proc.stdout.readlines():
2645 # return_code=proc.returncode
2647 # except subprocess.TimeoutExpired:
2649 # f=open(log_file,'a')
2650 # for line in proc.stdout.readlines():
2653 # f.write('Simulation stoppe because of TimeOut')
2660 def readObsnamesfile(__filenames=None):
2661 "Provisionnel pour lire les noms des observations d'un fichier"
2662 for __filename in __filenames:
2664 df = pandas.read_csv(__filename, sep = ";", header =0) #so that repeated names are not modified
2665 obsnames_infile = df.loc[0, :].values.tolist()[2:]
2667 df_tmp = pandas.read_csv(__filename, sep = ";", header =1)
2668 with open(__filename, 'r') as infile:
2669 readie=csv.reader(infile, delimiter=';')
2670 with open('tmp_obs_file.csv', 'wt', newline='') as output:
2671 outwriter=csv.writer(output, delimiter=';')
2672 list_row = ['#'] + df_tmp.columns
2673 outwriter.writerow(list_row)
2675 outwriter.writerow(row)
2677 df = pandas.read_csv('tmp_obs_file.csv', sep = ";", header =0) #so that repeated names are not modified
2678 df = df.iloc[1: , :]
2679 obsnames_infile = df.iloc[0, :].values.tolist()[2:]
2682 os.remove('tmp_obs_file.csv')
2685 return obsnames_infile
2687 def readDatesHours(__filename=None):
2688 "Provisionnel pour lire les heures et les dates d'une simulation dynamique"
2689 df = pandas.read_csv(__filename, sep = ";", header =1)
2690 dates = numpy.array(df['Date'])
2691 hours = numpy.array(df['Hour'])
2692 return (dates,hours)
2694 def readMultiDatesHours(__filenames=None):
2695 "Provisionnel pour lire les heures et les dates d'une simulation dynamique dans le cas multi-observations"
2698 for __filename in __filenames:
2699 df = pandas.read_csv(__filename, sep = ";", header =1)
2700 dates = numpy.array(df[df.columns[0]])
2701 hours = numpy.array(df[df.columns[1]])
2702 Multidates.append( dates )
2703 Multihours.append( hours )
2704 return (Multidates,Multihours)
2706 def find_nearest_index(array, value):
2707 "Trouver le point de la simulation qui se rapproche le plus aux observations"
2708 array = numpy.asarray(array)
2709 idx = (numpy.abs(array - value)).argmin()
2710 return [idx, array[idx]]
2711 # ==============================================================================
2712 if __name__ == "__main__":
2713 print('\n AUTODIAGNOSTIC\n ==============\n')
2714 print(" Module name..............: %s"%name)
2715 print(" Module version...........: %s"%version)
2716 print(" Module year..............: %s"%year)
2718 print(" Configuration attributes.:")
2719 for __k, __v in _configuration.items():
2720 print("%26s : %s"%(__k,__v))