COPASI API  4.16.103
SBMLImporter.h
Go to the documentation of this file.
1 // Copyright (C) 2010 - 2013 by Pedro Mendes, Virginia Tech Intellectual
2 // Properties, Inc., University of Heidelberg, and The University
3 // of Manchester.
4 // All rights reserved.
5 
6 // Copyright (C) 2008 - 2009 by Pedro Mendes, Virginia Tech Intellectual
7 // Properties, Inc., EML Research, gGmbH, University of Heidelberg,
8 // and The University of Manchester.
9 // All rights reserved.
10 
11 // Copyright (C) 2004 - 2007 by Pedro Mendes, Virginia Tech Intellectual
12 // Properties, Inc. and EML Research, gGmbH.
13 // All rights reserved.
14 
15 #ifndef SBMLIMPORTER_H__
16 #define SBMLIMPORTER_H__
17 
18 #include <string>
19 #include <map>
20 #include <set>
21 #include <utility>
22 #include "sbml/math/ASTNode.h"
23 
26 #include "copasi/model/CModel.h"
27 
28 class SBMLDocument;
29 class CCompartment;
30 class CMetab;
31 class CReaction;
32 class Reaction;
33 class Species;
34 class Model;
35 class Compartment;
36 class ConverterASTNode;
37 class Parameter;
38 class FunctionDefinition;
39 class SBase;
40 class CProcessReport;
41 class Rule;
42 class CListOfLayouts;
43 
45 {
46 protected:
47  static
48  C_FLOAT64 round(const C_FLOAT64 & x);
49 
50  static
51  bool areApproximatelyEqual(const double & x, const double & y, const double & t = 1e-9);
52 
53 protected:
54  std::set<unsigned int> mIgnoredSBMLMessages;
55  std::map<std::string, CMetab*> speciesMap;
63  unsigned int mLevel;
64  unsigned int mOriginalLevel;
65  unsigned int mVersion;
66  std::map<CFunction*, std::string> sbmlIdMap;
67  std::set<std::string> mUsedFunctions;
70  std::map<std::string, std::string> mFunctionNameMapping;
71  std::set<std::string> mDivisionByCompartmentReactions;
73  unsigned C_INT32 mImportStep;
74  size_t mhImportStep;
75  unsigned C_INT32 mTotalSteps;
76  std::map<Species*, Compartment*> mSubstanceOnlySpecies;
77  std::set<std::string> mFastReactions;
80  std::vector<std::string> mIgnoredParameterUnits;
81  std::map<const ASTNode*, CChemEqElement* > mStoichiometricExpressionMap;
83  std::set<const Parameter*> mPotentialAvogadroNumbers;
86  // this map maps a delay expression to the global parameter id it has been
87  // replaced with
88  std::map<std::string, std::string> mDelayNodeMap;
89  std::set<std::string> mUsedSBMLIds;
92  std::map<std::string, std::string> mKnownCustomUserDefinedFunctions;
93  std::map<std::string, std::string> mKnownInitalValues;
94 
95 #if LIBSBML_VERSION >= 40100
96  // this map is used for storing the parameters that are used as factors that have to be applied to the multiplicities
97  // of the chemical equation elements
98  const CModelValue* mpModelConversionFactor;
99  // we only store the id of the species as the value and use this value as the key into the mSpeciesConversionParameterMap below
100  std::map<CChemEqElement*, std::pair<std::string, CChemEq::MetaboliteRole> > mChemEqElementSpeciesIdMap;
101  // yet another map that stores conversion parameters per species id
102  // This will speed up the assignment of
103  std::map<std::string, const CModelValue*> mSpeciesConversionParameterMap;
104  // and yet another map that maps SBML parameter keys to COPASI model values
105  std::map<std::string, const CModelValue*> mSBMLIdModelValueMap;
106 
107  // in this set we store the ids of all SBML species references
108  // so that we can later check if a reference to a species reference is used
109  // in a mathematical expression
110  std::map<std::string, double> mSBMLSpeciesReferenceIds;
111 
112  bool mRateRuleForSpeciesReferenceIgnored;
113  bool mEventAssignmentForSpeciesReferenceIgnored;
114  bool mConversionFactorFound;
115 #if LIBSBML_VERSION >= 40200
116  bool mEventPrioritiesIgnored;
117  bool mInitialTriggerValues;
118  bool mNonPersistentTriggerFound;
119 #endif // LIBSBML_VERSION >= 40200
120 
121 #endif // LIBSBML_VERSION >= 40100
122 
123  /**
124  * Creates and returns a COPASI CModel from the SBMLDocument given as argument.
125  */
126  CModel* createCModelFromSBMLDocument(SBMLDocument * doc,
127  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
128 
129  /**
130  * Creates and returns a COPASI CFunction from the SBML FunctionDefinition
131  * given as argument.
132  */
133  CFunction* createCFunctionFromFunctionDefinition(const FunctionDefinition* sbmlFunction,
134  CFunctionDB* pTmpFunctionDB,
135  Model* pSBMLModel,
136  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
137 
138  CFunction* createCFunctionFromFunctionTree(const FunctionDefinition* pSBMLFunction,
139  Model* pSBMLModel,
140  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
141 
142  /**
143  * Creates and returns a COPASI CCompartment from the SBML Compartment
144  * given as argument.
145  */
146  CCompartment* createCCompartmentFromCompartment(const Compartment* sbmlComp,
147  CModel* copasiModel,
148  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
149  const Model* pSBMLModel);
150 
151  /**
152  * Creates and returns a COPASI CMetab from the given SBML Species object.
153  */
154  CMetab* createCMetabFromSpecies(const Species* sbmlSpecies,
155  CModel* copasiModel,
156  CCompartment* copasiCompartment,
157  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
158  const Model* pSBMLModel);
159 
160  /**
161  * Checks if no id is used in more than one Assignment and RateRule.
162  */
163  void areRulesUnique(const Model* copasiMode);
164 
165  /**
166  * Imports the given Rule if COPASI supports this kind of Rule, otherwise a warning is created.
167  */
168  void importSBMLRule(const Rule* sbmlRule,
169  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
170  Model* pSBMLModel);
171 
172  /**
173  * Imports the given AssignmentRule which is for a global parameter.
174  */
175  void importRuleForModelEntity(const Rule* rule,
176  CModelEntity* pMV,
177  CModelEntity::Status ruleType,
178  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
179  Model* pSBMLModel);
180 
181  /**
182  * Imports all events
183  */
184  void importEvents(Model* pSBMLModel,
185  CModel* pCopasiModel,
186  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
187 
188  /**
189  * Imports the given event.
190  */
191  void importEvent(const Event* pEvent,
192  Model* pSBMLModel,
193  CModel* pCopasiModel,
194  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
195 
196  /**
197  * Imports the given RateRule if COPASI supports this kind of RateRule, otherwise a warning is created.
198  */
199  void importRule(const Rule* rule,
200  CModelEntity::Status ruleType,
201  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
202  Model* pSBMLModel);
203 
204  /**
205  * Recurses an ASTNode tree and gets all SBML Ids in the tree.
206  * The ids are stored in the given set.
207  */
208  void getIdsFromNode(const ASTNode* pNode,
209  std::set<std::string>& idSet);
210 
211  /**
212  * Checks the expression for a give rate or assignment rule for
213  * consistency. This basically means it checks that no id present in the
214  * expression is the target for one of the following rate or assignment
215  * rules.
216  */
217  void checkRuleMathConsistency(const Rule* pRule,
218  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
219 
220  /**
221  * Creates and returns a COPASI CModelValue from the given SBML Parameter object.
222  */
223  CModelValue* createCModelValueFromParameter(const Parameter* sbmlParameter,
224  CModel* copasiModel,
225  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
226 
227  /**
228  * Creates and returns a COPASI CReaction object from the given SBML
229  * Reaction object.
230  */
231  CReaction* createCReactionFromReaction(Reaction* sbmlReaction,
232  Model* sbmlModel,
233  CModel* cmodel, std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
234  CFunctionDB* pTmpFunctionDB);
235 
236  /**
237  * Creates a map of each parameter of the function definition and its
238  * corresponding parameter in the function call.
239  */
240  std::map<std::string , ASTNode*> createBVarMap(const ASTNode* uDefFunction,
241  const ASTNode* function);
242 
243  /**
244  * Returns the user defined SBML function definition that belongs to the given
245  * name, or NULL if none can be found.
246  */
247  const FunctionDefinition* getFunctionDefinitionForName(const std::string name,
248  const Model* model);
249 
250  /**
251  * Replaces the variables in a function definition with the actual function
252  * parameters that were used when the function was called. The function returns
253  * a pointer to the ConverterAST node with the replaced variables.
254  */
255  // ConverterASTNode* replaceBvars(const ASTNode* node,
256  // std::map<std::string, ASTNode*> bvarMap);
257 
258  /**
259  * This function replaces the AST_FUNCTION_POWER ASTNodes in a ASTNode tree
260  * with the AST_POWER node.
261  */
262  void replacePowerFunctionNodes(ASTNode* node);
263 
264  /**
265  * This functions replaces all species nodes for species that are in the substanceOnlySpeciesVector.
266  * With the node multiplied by the volume of the species compartment.
267  void replaceSubstanceOnlySpeciesNodes(ConverterASTNode* node, const std::map<Species*, Compartment*>& substanceOnlySpecies);
268  */
269 
270  /**
271  * Returns the copasi LengthUnit corresponding to the given SBML length
272  * UnitDefinition.
273  */
274  std::pair<CModel::LengthUnit, bool> handleLengthUnit(const UnitDefinition* uDef);
275 
276  /**
277  * Returns the copasi AreaUnit corresponding to the given SBML area
278  * UnitDefinition.
279  */
280  std::pair<CModel::AreaUnit, bool> handleAreaUnit(const UnitDefinition* uDef);
281 
282  /**
283  * Returns the copasi VolumeUnit corresponding to the given SBML Volume
284  * UnitDefinition.
285  */
286  std::pair<CModel::VolumeUnit, bool> handleVolumeUnit(const UnitDefinition* uDef);
287 
288  /**
289  * Returns the COPASI QuantityUnit corresponding to the given SBML
290  * Substance UnitDefinition.
291  */
292  std::pair<CModel::QuantityUnit, bool> handleSubstanceUnit(const UnitDefinition* uDef);
293 
294  /**
295  * Returns the COPASI TimeUnit corresponding to the given SBML Time
296  * UnitDefinition.
297  */
298  std::pair<CModel::TimeUnit, bool> handleTimeUnit(const UnitDefinition* uDef);
299 
300  /**
301  * Replaces all occurrences of the log function with two arguments by
302  * a division of two separate calls to log.
303  */
304  void replaceLog(ConverterASTNode* sourceNode);
305 
306  /**
307  * Replaces all occurrences of the root function with two arguments by
308  * a call to the power function with the inverse of the first argument.
309  */
310  void replaceRoot(ConverterASTNode* sourceNode);
311 
312  /**
313  * Replaces the ids of named nodes in an ASTNode tree with
314  * the correspondingCopasi Common Names.
315  */
316  bool sbmlId2CopasiCN(ASTNode* pNode,
317  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
318  CCopasiParameterGroup& pParamGroup);
319 
320  /**
321  * Upon import a function object might change its name due to naming conflicts in the function
322  * database. So each ASTNode tree has its call node names replaced
323  * before it will be processed further.
324  */
325  void replaceCallNodeNames(ASTNode* pNode);
326 
327  /**
328  * This function replaces calls to the delay function in an ASTNode tree
329  * by a node that references a new global parameter which the function
330  * creates. The global parameter gets an expression which corresponds to the
331  * delay call.
332  * This is necessary because all knetic laws in COPASI are function calls and
333  * function definitions should not contain a call to delay.
334  */
336  Model* pModel,
337  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
338  Reaction* pSBMLReaction,
339  std::map<std::string, std::string>& localReplacementMap
340  );
341 
342  /**
343  * If we replace delay nodes within kineitc laws, we have to make sure that
344  * there is no reference to a local parameter within the replaced
345  * delay node because that would mean that we end up with a reference to a
346  * local parameter in the rule for the delay replacement which is not allowed
347  * in SBML.
348  * Therefore we have to convert all local parameters which occur within a
349  * delay call into global parameters.
350  * This method finds all local parameters that have to be converted to global
351  * parameters and it already creates the necessary global parameters.
352  * The localReplacementMap returns a mapping between the is of the original
353  * parameter and the id of the new parameter it will be replaced with.
354  * This map is used in a second step to actually replace the nodes in the
355  * expression.
356  */
357  void find_local_parameters_in_delay(ASTNode* pNode,
358  Reaction* pSBMLReaction,
359  Model* pModel,
360  std::map<std::string, std::string>& localReplacementMap,
361  const std::set<std::string>& localIds,
362  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap
363  );
364 
365  /**
366  * This method gets an ASTNode and a map between old node names and new node
367  * names. All AST_NAME nodes with an "old" name are replaced by a node with
368  * the "new" name.
369  */
370  void replace_name_nodes(ASTNode* pNode, const std::map<std::string, std::string>& replacementMap);
371 
372  /**
373  * The data for a CEvaluationNodeObject needs to have the common
374  * name of the model it refers to as its data. Since this model is
375  * only known via a pointer in the SBMLImporter at the time of
376  * import, all AST_NAME_TIME nodes that are imported need to have
377  * their name replaced by the common name of this model.
378  * Starting with SBML Level 3 this also applies to the avogadro number.
379  */
380  void replaceTimeAndAvogadroNodeNames(ASTNode* pNode);
381 
382  /**
383  * COPASI can not handle the delay function yet, so if it is used in some expression, we
384  * have to abort reading the file.
385  */
387 
388  /**
389  * In a preprocessing step each expression tree that is imported,
390  * e.g. for function definitions, rules, event or kinetic laws
391  * is preprocessed to replace some of the nodes data.
392  * See also replaceCallNodeNames and replaceTimeNodeNames.
393  */
394  void preprocessNode(ConverterASTNode* pNode,
395  Model* pSBMLModel,
396  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
397  Reaction* pSBMLReaction = NULL);
398 
399  CFunction* findCorrespondingFunction(const CExpression * pExpression, const CReaction* reaction);
400 
401 public:
402  static bool areEqualFunctions(const CFunction* pFun, const CFunction* pFun2);
403 
404  /**
405  * Compares to CEvaluationNode based subtrees recursively.
406  */
407  static bool areEqualSubtrees(const CEvaluationNode* pNode1, const CEvaluationNode* pNode2);
408 
409 protected:
410  std::vector<CEvaluationNodeObject*>* isMassAction(const CEvaluationTree* pTree,
411  const CChemEq& chemicalEquation,
412  const CEvaluationNodeCall* pCallNode = NULL);
413 
414  /**
415  * Checks if the given node is an object node that represents a parameter
416  * or a model value or a function which has a single parameter and a single node which also represents a parameter.
417  */
419  CModel* pModel,
420  CFunctionDB* pFunctionDB);
421 
422  std::vector<CEvaluationNodeObject*>* isMassActionExpression(const CEvaluationNode* pRootNode,
423  const CChemEq& chemicalEquation);
424 
425  std::vector<CEvaluationNodeObject*>* isMassActionFunction(const CFunction* pFun,
426  const CChemEq& chemicalEquation,
427  const std::vector<std::vector< std::string > >& functionArgumentCNs);
428 
429  /**
430  * This function takes a node and tries to find out whether the tree under this node consists
431  * only of multiply operators and object nodes.
432  * The arguments to the multiply operators are returned.
433  */
434  void separateProductArguments(const CEvaluationNode* pRootNode,
435  std::vector<const CEvaluationNode*>& arguments);
436 
437  void doMapping(CReaction* pCopasiReaction, const CEvaluationNodeCall* pCallNode);
438 
439  bool isSimpleFunctionCall(const CEvaluationNode* pRootNode);
440 
441  void setCorrectUsage(CReaction* pCopasiReaction, const CEvaluationNodeCall* pCallNode);
442 
443  /**
444  * Checks if a given tree is multiplied by a compartment identifier.
445  * If so, a copy of the tree is returned in which the multiplication has been removed.
446  */
447  ConverterASTNode* isMultipliedByVolume(const ASTNode* node, const std::string& compartmentSBMLId);
448 
450  const std::map<std::string, std::string>& replacementMap);
451 
453  const std::vector<std::vector<std::string > >& functionArgumentCNs);
454 
456 
457  /**
458  * Checks if the volume with the given CN is one of the parameters to the given call node.
459  */
460  bool containsVolume(const ASTNode* pNode, const std::string& compartmentCN);
461 
462  /**
463  * Finds all functions that are used and removes those that are not.
464  */
465  bool removeUnusedFunctions(CFunctionDB* pTmpFunctionDB,
466  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
467 
468  /**
469  * Finds all functions calls directly or indirectly used in a function
470  * tree.
471  */
472  void findFunctionCalls(const CEvaluationNode* pNode, std::set<std::string>& functionNameSet);
473 
474  /**
475  * Heuristically checks whether a model was meant to be simulated stochastically.
476  * If the substance units are set to items and all reaction are irreversible this function
477  * will return true;
478  */
479  bool isStochasticModel(const Model* pSBMLModel);
480 
481  /**
482  * If an initial expression uses time, we have to import it as initial time
483  * instead.
484  * This method takes an AST node and converts all time nodes to object nodes
485  * that have the common name of the time as the name.
486  */
487  void replace_time_with_initial_time(ASTNode* pNode, const CModel* pCOPASIModel);
488 
489  void replaceObjectNames(ASTNode* pNode,
490  const std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
491  bool initialExpression = false);
492 
493  /**
494  * For function definitions that use the time symbol we have to make this a
495  * variable that is passed to the function instead.
496  * The function recursively goes through the AST tree rooted in root and
497  * changes all time nodes to variable nodes with name newNodeName.
498  * If a time node has been found, the function return true, otherwise false
499  * is returned.
500  */
501  bool replaceTimeNodesInFunctionDefinition(ASTNode* root, std::string newNodeName);
502 
503  /**
504  * This function replaces function calls to all functions listed in mExplicitelyTimeDependentFunctionDefinitions
505  * with the same call but an additional parameter which is the time.
506  * This replacement includes all model entities that have a mathematical
507  * expression, so depending on the version of SBML, this would include:
508  * initial assignments, rules, constraints, kinetic laws and events.
509  * The corresponding replacement for the function definitions is done in
510  * replaceTimeNodesInFunctionDefinition.
511  */
512  void replaceTimeDependentFunctionCalls(ASTNode* root);
513 
514  /**
515  * Sets the initial values on compartments, metabolites and model values if
516  * those initial values have been set in the SBML model.
517  * Otherwise the routine checks if a rule or an initial assignment has been set for the entity.
518  * If the entity has not been set in any way, an error message is created.
519  */
520  bool setInitialValues(CModel* pModel,
521  const std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
522 
523  void checkElementUnits(const Model* pSBMLModel,
524  CModel* pCopasiModel,
525  int level,
526  int version);
527 
528  /**
529  * If the given UnitDefinition can be converted to a form of litre, the
530  * function return the UnitDefinition in litre, otherwise NULL is returned.
531  */
532  static Unit* convertSBMLCubicmetresToLitres(const Unit* pU);
533 
534  /**
535  * This function normalizes the multiplier to be within the range 1.0 <=
536  * multiplier < 10.0.
537  */
538  static void normalizeSBMLUnit(Unit* pU);
539 
540  /**
541  * Imports all initial assignments if there are any.
542  */
543  void importInitialAssignments(Model* pSBMLModel, std::map<CCopasiObject*, SBase*>& copasi2sbmlMap, const CModel* pCOPASIModel);
544 
545  /**
546  * This method evaluates all stoichiometric expressions and sets them as
547  * constants on the CChemEqElement.
548  */
549  void applyStoichiometricExpressions(std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
550  Model* pSBMLModel);
551 
552  /**
553  * Creates a function definition for the delay function.
554  */
556 
557  /**
558  * This method goes through the list of global parameters and tries to find
559  * a parameter that could correspond to avogadros number.
560  */
561  void findAvogadroConstant(Model* pSBMLModel,
562  double factor);
563 
564  /**
565  * This method replaces references to the id of species which have the
566  * hasOnlySubstanceUnits flag set with the reference divided by avogadros
567  * number.
568  * The method tries to determine if there already is a multiplication with
569  * avogadros number and removes this multiplication rather than adding a new division.
570  */
572  Model* pSBMLModel,
573  double factor,
574  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
575 
576  /**
577  * This method creates a global parameter the represents the factor that is
578  * used to convert a particle number into the amount units set on the
579  * model.
580  * The parameter is only created if it is needed and after exporting the
581  * model, the parameter is deleted from the COPASI model again.
582  */
583  void createHasOnlySubstanceUnitFactor(Model* pSBMLModel,
584  double factor,
585  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
586 
587  /**
588  * Multiplies all species nodes that belong to species with the
589  * hasSubstanceOnlyUnits flag set with the volume of the compartment that
590  * the species belongs to.
591  * This is only done for kineticLaw, all other mathematical expressions
592  * import those references as particle number nodes divided by the
593  * quantity2unit factor.
594  */
596 
597  /**
598  * Imports the MIRIAM annotation from the given SBML object and sets it on
599  * the given COPASI object.
600  */
601  bool importMIRIAM(const SBase* pSBMLObject, CCopasiObject* pCOPASIObject);
602 
603  /**
604  * Starting with SBML Level 2 Version 4, function definitions no longer need to
605  * be ordered, i.e. a function definition may refer to another function
606  * definition that is defined somewhere further down in the file.
607  * So we have to import the function definitions in the correct order.
608  */
609  CFunctionDB* importFunctionDefinitions(Model* pSBMLModel,
610  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
611 
612  /**
613  * static method that finds all direct function dependencies of a given
614  * function definition.
615  */
616  static void findDirectDependencies(const FunctionDefinition* pFunDef,
617  std::map<const FunctionDefinition*, std::set<std::string> >& dependencies);
618 
619  /**
620  * static method that recursively finds all direct function dependencies of the
621  * expression rooted at the given node.
622  */
623  static void findDirectDependencies(const ASTNode* pNode,
624  std::set<std::string>& dependencies);
625 
626 public:
627  SBMLImporter();
628  ~SBMLImporter();
629 
630  CModel* readSBML(std::string filename,
631  CFunctionDB* funDB,
632  SBMLDocument*& pSBMLDocument,
633  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
634  CListOfLayouts *& prLol,
635  CCopasiDataModel* pDataModel);
636 
637  //CModel* readSBML(std::string filename, CFunctionDB* funDB, SBMLDocument *& pSBMLDocument, std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
638 
639  CModel* parseSBML(const std::string& sbmlDocumentText,
640  CFunctionDB* funDB,
641  SBMLDocument *& pSBMLDocument,
642  std::map<CCopasiObject*, SBase*>& copasi2sbmlmap,
643  CListOfLayouts *& prLol,
644  CCopasiDataModel* pDataModel);
645 
646  //CModel* parseSBML(const std::string& sbmlDocumentText, CFunctionDB* funDB, SBMLDocument*& pSBMLDocument, std::map<CCopasiObject*, SBase*>& copasi2sbmlmap);
647 #ifdef COPASI_DEBUG
648  static void printMap(const std::map<CCopasiObject*, SBase*>& map);
649 #endif // COPASI_DEBUG
650 
651  void restoreFunctionDB();
652 
653  /**
654  * This call deletes an existing COPASI model.
655  * The method can e.g. be called to clean up if an import fails.
656  */
657  void deleteCopasiModel();
658 
659  void setImportHandler(CProcessReport* pHandler);
660 
661  /**
662  * Enhanced method to identify identical SBML unit definitions.
663  * This method uses the areIdentical method from libSBML, but if the method
664  * return false, it does some extra checks.
665  * Right now it check for example if two volumes, one given in litre and one
666  * given in cubic meters are identical.
667  */
668  static bool areSBMLUnitDefinitionsIdentical(const UnitDefinition* pUdef1, const UnitDefinition* pUdef2);
669 
670  /**
671  * This method takes the id of a unit as it can appear in an SBML file, and
672  * returns a new UnitDefinition object for that id.
673  */
674  static UnitDefinition* getSBMLUnitDefinitionForId(const std::string& unitId,
675  const Model* pSBMLModel);
676 
678 
679  /**
680  * Returns the flag that determines whether COPASI MIRIAM annotation is
681  * imported if it is present.
682  */
683  bool getImportCOPASIMIRIAM() const;
684 
685  /**
686  * Sets the flag that determines whether COPASI MIRIAM annotation is
687  * imported if it is present.
688  */
689  void setImportCOPASIMIRIAM(bool import);
690 
691  /**
692  * This method takes an AST node and a set of ids and returns the first id
693  * from the set it finds in the AST tree.
694  * This is e.g. used to check if expression in L2V1 contain references to reaction ids.
695  */
696  std::string findIdInASTTree(const ASTNode* pMath, const std::set<std::string>& reactionIds);
697  std::string findIdInASTTree(const ASTNode* pMath, const std::map<std::string, double>& reactionIds);
698 
699  /**
700  * This method divides the given expression by the given object and returns a new expression.
701  * The caller is responsible for freeing the memory for the new expression.
702  */
703  static CEvaluationNode* divideByObject(const CEvaluationNode* pOrigNode, const CCopasiObject* pObject);
704 
705  /**
706  * This method reads the notes from an arbitrate SBase object
707  * and set them on the given CAnnotation instance.
708  */
709  static bool importNotes(CAnnotation* pAnno, const SBase* pSBase);
710 
711 #if LIBSBML_VERSION >= 40100
712  /**
713  * This method check if a unit has been set on a number node.
714  * If such a node is found in the tree, true is returned.
715  */
716  static bool checkForUnitsOnNumbers(const ASTNode* pNode);
717 
718  /**
719  * This method checks if there are conversion factors that need to be applied to
720  * ChemicalEquationElements and applies them.
721  */
722  void applyConversionFactors();
723 
724  /**
725  * Goes through all SBML reactions and collects the ids of all species references.
726  */
727  static void updateSBMLSpeciesReferenceIds(const Model* pModel, std::map<std::string, double>& ids);
728 
729 #endif // LIBSBML_VERSION
730 };
731 
732 #endif // SBMLIMPORTER_H__
void multiplySubstanceOnlySpeciesByVolume(ConverterASTNode *pNode)
void setImportCOPASIMIRIAM(bool import)
void createDelayFunctionDefinition()
std::pair< CModel::TimeUnit, bool > handleTimeUnit(const UnitDefinition *uDef)
std::set< unsigned int > mIgnoredSBMLMessages
Definition: SBMLImporter.h:54
size_t mhImportStep
Definition: SBMLImporter.h:74
void replace_delay_nodes(ConverterASTNode *pNode, Model *pModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, Reaction *pSBMLReaction, std::map< std::string, std::string > &localReplacementMap)
static void normalizeSBMLUnit(Unit *pU)
void restoreFunctionDB()
void importSBMLRule(const Rule *sbmlRule, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, Model *pSBMLModel)
static C_FLOAT64 round(const C_FLOAT64 &x)
bool mAvogadroCreated
Definition: SBMLImporter.h:84
void checkRuleMathConsistency(const Rule *pRule, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
bool mUnitOnNumberFound
Definition: SBMLImporter.h:61
void importRule(const Rule *rule, CModelEntity::Status ruleType, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, Model *pSBMLModel)
bool mUnsupportedRateRuleFound
Definition: SBMLImporter.h:59
unsigned int mLevel
Definition: SBMLImporter.h:63
bool containsVolume(const ASTNode *pNode, const std::string &compartmentCN)
std::map< std::string, std::string > mKnownInitalValues
Definition: SBMLImporter.h:93
CModel * mpCopasiModel
Definition: SBMLImporter.h:69
void replaceTimeDependentFunctionCalls(ASTNode *root)
void find_local_parameters_in_delay(ASTNode *pNode, Reaction *pSBMLReaction, Model *pModel, std::map< std::string, std::string > &localReplacementMap, const std::set< std::string > &localIds, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
std::map< std::string, std::string > mFunctionNameMapping
Definition: SBMLImporter.h:70
std::set< std::string > mReactionsWithReplacedLocalParameters
Definition: SBMLImporter.h:78
std::vector< CEvaluationNodeObject * > * isMassActionFunction(const CFunction *pFun, const CChemEq &chemicalEquation, const std::vector< std::vector< std::string > > &functionArgumentCNs)
void importInitialAssignments(Model *pSBMLModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlMap, const CModel *pCOPASIModel)
void findAvogadroConstant(Model *pSBMLModel, double factor)
void replace_name_nodes(ASTNode *pNode, const std::map< std::string, std::string > &replacementMap)
static Unit * convertSBMLCubicmetresToLitres(const Unit *pU)
std::map< Species *, Compartment * > mSubstanceOnlySpecies
Definition: SBMLImporter.h:76
void deleteCopasiModel()
void setImportHandler(CProcessReport *pHandler)
std::map< std::string, ASTNode * > createBVarMap(const ASTNode *uDefFunction, const ASTNode *function)
CProcessReport * mpImportHandler
Definition: SBMLImporter.h:72
CFunction * createCFunctionFromFunctionDefinition(const FunctionDefinition *sbmlFunction, CFunctionDB *pTmpFunctionDB, Model *pSBMLModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
void renameMassActionParameters(CEvaluationNodeCall *pCallNode)
bool isDelayFunctionUsed(ConverterASTNode *pNode)
bool mDelayFound
Definition: SBMLImporter.h:82
void getIdsFromNode(const ASTNode *pNode, std::set< std::string > &idSet)
CModel * createCModelFromSBMLDocument(SBMLDocument *doc, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
void replaceObjectNames(ASTNode *pNode, const std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, bool initialExpression=false)
bool mUnsupportedRuleFound
Definition: SBMLImporter.h:58
void replace_time_with_initial_time(ASTNode *pNode, const CModel *pCOPASIModel)
CFunctionDB * importFunctionDefinitions(Model *pSBMLModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
unsigned int mVersion
Definition: SBMLImporter.h:65
void importEvent(const Event *pEvent, Model *pSBMLModel, CModel *pCopasiModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
#define C_INT32
Definition: copasi.h:90
std::set< const Parameter * > mPotentialAvogadroNumbers
Definition: SBMLImporter.h:83
void replaceLog(ConverterASTNode *sourceNode)
Definition: CMetab.h:178
void replaceRoot(ConverterASTNode *sourceNode)
CMetab * createCMetabFromSpecies(const Species *sbmlSpecies, CModel *copasiModel, CCompartment *copasiCompartment, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, const Model *pSBMLModel)
CCopasiObject * isConstantFlux(const CEvaluationNode *pRoot, CModel *pModel, CFunctionDB *pFunctionDB)
void checkElementUnits(const Model *pSBMLModel, CModel *pCopasiModel, int level, int version)
void createHasOnlySubstanceUnitFactor(Model *pSBMLModel, double factor, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
std::set< std::string > mUsedFunctions
Definition: SBMLImporter.h:67
std::set< std::string > mUsedSBMLIds
Definition: SBMLImporter.h:89
CProcessReport * getImportHandlerAddr()
CFunction * createCFunctionFromFunctionTree(const FunctionDefinition *pSBMLFunction, Model *pSBMLModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
void separateProductArguments(const CEvaluationNode *pRootNode, std::vector< const CEvaluationNode * > &arguments)
ConverterASTNode * isMultipliedByVolume(const ASTNode *node, const std::string &compartmentSBMLId)
unsigned C_INT32 mImportStep
Definition: SBMLImporter.h:73
std::pair< CModel::QuantityUnit, bool > handleSubstanceUnit(const UnitDefinition *uDef)
static CEvaluationNode * divideByObject(const CEvaluationNode *pOrigNode, const CCopasiObject *pObject)
void setCorrectUsage(CReaction *pCopasiReaction, const CEvaluationNodeCall *pCallNode)
std::vector< CEvaluationNodeObject * > * isMassAction(const CEvaluationTree *pTree, const CChemEq &chemicalEquation, const CEvaluationNodeCall *pCallNode=NULL)
std::map< std::string, std::string > mKnownCustomUserDefinedFunctions
Definition: SBMLImporter.h:92
CCopasiDataModel * mpDataModel
Definition: SBMLImporter.h:68
std::vector< CEvaluationNodeObject * > * isMassActionExpression(const CEvaluationNode *pRootNode, const CChemEq &chemicalEquation)
bool mAvogadroSet
Definition: SBMLImporter.h:91
CFunctionDB * functionDB
Definition: SBMLImporter.h:56
static bool areEqualSubtrees(const CEvaluationNode *pNode1, const CEvaluationNode *pNode2)
static bool areEqualFunctions(const CFunction *pFun, const CFunction *pFun2)
bool setInitialValues(CModel *pModel, const std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
CModel * parseSBML(const std::string &sbmlDocumentText, CFunctionDB *funDB, SBMLDocument *&pSBMLDocument, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, CListOfLayouts *&prLol, CCopasiDataModel *pDataModel)
bool mIncompleteModel
Definition: SBMLImporter.h:57
CFunction * findCorrespondingFunction(const CExpression *pExpression, const CReaction *reaction)
bool mAssignmentToSpeciesReferenceFound
Definition: SBMLImporter.h:62
void areRulesUnique(const Model *copasiMode)
CEvaluationNode * variables2objects(const CEvaluationNode *pOrigNode, const std::map< std::string, std::string > &replacementMap)
void preprocessNode(ConverterASTNode *pNode, Model *pSBMLModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, Reaction *pSBMLReaction=NULL)
bool sbmlId2CopasiCN(ASTNode *pNode, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, CCopasiParameterGroup &pParamGroup)
bool isSimpleFunctionCall(const CEvaluationNode *pRootNode)
std::map< const ASTNode *, CChemEqElement * > mStoichiometricExpressionMap
Definition: SBMLImporter.h:81
std::set< std::string > mDivisionByCompartmentReactions
Definition: SBMLImporter.h:71
bool removeUnusedFunctions(CFunctionDB *pTmpFunctionDB, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
CEvaluationTree * createExpressionFromFunction(const CFunction *pFun, const std::vector< std::vector< std::string > > &functionArgumentCNs)
bool isStochasticModel(const Model *pSBMLModel)
void findFunctionCalls(const CEvaluationNode *pNode, std::set< std::string > &functionNameSet)
std::pair< CModel::AreaUnit, bool > handleAreaUnit(const UnitDefinition *uDef)
static UnitDefinition * getSBMLUnitDefinitionForId(const std::string &unitId, const Model *pSBMLModel)
std::set< std::string > mFastReactions
Definition: SBMLImporter.h:77
void doMapping(CReaction *pCopasiReaction, const CEvaluationNodeCall *pCallNode)
std::map< std::string, std::string > mDelayNodeMap
Definition: SBMLImporter.h:88
void replaceCallNodeNames(ASTNode *pNode)
CCompartment * createCCompartmentFromCompartment(const Compartment *sbmlComp, CModel *copasiModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, const Model *pSBMLModel)
std::set< std::string > mExplicitelyTimeDependentFunctionDefinitions
Definition: SBMLImporter.h:79
std::string findIdInASTTree(const ASTNode *pMath, const std::set< std::string > &reactionIds)
#define C_FLOAT64
Definition: copasi.h:92
const FunctionDefinition * getFunctionDefinitionForName(const std::string name, const Model *model)
bool mImportCOPASIMIRIAM
Definition: SBMLImporter.h:85
static bool areApproximatelyEqual(const double &x, const double &y, const double &t=1e-9)
bool replaceTimeNodesInFunctionDefinition(ASTNode *root, std::string newNodeName)
void importEvents(Model *pSBMLModel, CModel *pCopasiModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
The class for handling a chemical kinetic function.
Definition: CFunction.h:29
std::map< CFunction *, std::string > sbmlIdMap
Definition: SBMLImporter.h:66
static void findDirectDependencies(const FunctionDefinition *pFunDef, std::map< const FunctionDefinition *, std::set< std::string > > &dependencies)
std::pair< CModel::LengthUnit, bool > handleLengthUnit(const UnitDefinition *uDef)
bool importMIRIAM(const SBase *pSBMLObject, CCopasiObject *pCOPASIObject)
CReaction * createCReactionFromReaction(Reaction *sbmlReaction, Model *sbmlModel, CModel *cmodel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, CFunctionDB *pTmpFunctionDB)
CModelValue * createCModelValueFromParameter(const Parameter *sbmlParameter, CModel *copasiModel, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
unsigned C_INT32 mTotalSteps
Definition: SBMLImporter.h:75
bool getImportCOPASIMIRIAM() const
void replacePowerFunctionNodes(ASTNode *node)
std::map< std::string, CMetab * > speciesMap
Definition: SBMLImporter.h:55
std::vector< std::string > mIgnoredParameterUnits
Definition: SBMLImporter.h:80
Definition: CModel.h:50
void replaceTimeAndAvogadroNodeNames(ASTNode *pNode)
CModel * readSBML(std::string filename, CFunctionDB *funDB, SBMLDocument *&pSBMLDocument, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, CListOfLayouts *&prLol, CCopasiDataModel *pDataModel)
static bool importNotes(CAnnotation *pAnno, const SBase *pSBase)
void replaceAmountReferences(ConverterASTNode *pNode, Model *pSBMLModel, double factor, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap)
std::pair< CModel::VolumeUnit, bool > handleVolumeUnit(const UnitDefinition *uDef)
void importRuleForModelEntity(const Rule *rule, CModelEntity *pMV, CModelEntity::Status ruleType, std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, Model *pSBMLModel)
unsigned int mOriginalLevel
Definition: SBMLImporter.h:64
void applyStoichiometricExpressions(std::map< CCopasiObject *, SBase * > &copasi2sbmlmap, Model *pSBMLModel)
static bool areSBMLUnitDefinitionsIdentical(const UnitDefinition *pUdef1, const UnitDefinition *pUdef2)
bool mUnsupportedAssignmentRuleFound
Definition: SBMLImporter.h:60
bool mUsedSBMLIdsPopulated
Definition: SBMLImporter.h:90