DianaICSolver: HafkeReactor.py

import diana 
try: 
    dmain=diana.GetDianaMain()
    # get Diana model manger 
    mmanger=dmain.GetModelManager();

    # create model 
    model=mmanger.CreateModel(diana.CAPE_CONTINUOUS, "HafkeReactor.so");
    model.Initialize(); 
    eso=model.GetActiveESO(); 
    esopar=eso.GetParameters(); 
    ri=dmain.CreateReportingInterface("basic");
    ri.SetComponentName("data"); 
    sfactory=dmain.GetSolverFactory(); 

    # create consistent initial conditions solver
    icsolver=sfactory.CreateSolver(diana.CAPE_IC, model, "icsolver.so");
    icsolver.Initialize();
    icpar=icsolver.GetParameters();
    icpar["VerboseLevel"].SetValue(0);
    icsolver.Solve();

    # create solver 
    solver=sfactory.CreateSolver(diana.CAPE_DAE, model, "ida.so"); 
    solver.Initialize(); 
    solpar=solver.GetParameters(); 
    solver.SetReportingInterface(ri); 
    ri.Add(solpar["T"]); 
    ri.Add(eso.GetStateVariables().ItemByName("c_[1]")); 
    ri.Add(eso.GetStateVariables().ItemByName("c_[2]")); 
    ri.Add(eso.GetStateVariables().ItemByName("c_[3]")); 
    ri.Add(eso.GetStateVariables().ItemByName("c_[4]")); 
    ri.Add(eso.GetStateVariables().ItemByName("temp")); 
    ri.Add(eso.GetStateVariables().ItemByName("tk")); 

    solpar["VerboseLevel"].SetValue(2);
    solpar["Tend"].SetValue(16000); 
    solver.Solve(); 

    # create directory and store results 
    import os
    if not os.path.isdir('results'): os.mkdir('results'); 
    ri.WriteData("results/HafkeReactor"); 
    ri.WriteDataMatlab("results/HafkeReactor.m"); 
except diana.ECapeUser, exc: 
    print exc