//-------------------------------------------------------------------
// $Id: m_param.cpp 1406 2009-08-21 08:37:56Z gems $
//
// Implementation of config and basic methods of TProfile class
//
// Rewritten from C to C++ by S.Dmytriyeva
// Copyright (C) 1995-2020 S.Dmytriyeva, D.Kulik
//
// This file is part of a GEM-Selektor library for thermodynamic
// modelling by Gibbs energy minimization
// Uses: GEM-Selektor GUI GUI DBMS library, gems/lib/gemvizor.lib
//
// This file may be distributed under the GPL v.3 license

//
// See http://gems.web.psi.ch/ for more information
// E-mail: gems2.support@psi.ch
//-------------------------------------------------------------------
//

#include <cstdint>
#include "m_unspace.h"
#include "m_gtdemo.h"
#include "m_syseq.h"
#include "m_proces.h"
#include "m_dualth.h"
#include "visor.h"
#include "nodearray_gui.h"
#include "GEMS3K/gdatastream.h"
#include "gemsgui_version.h"
#include "GEMS3K/gems3k_version.h"

#ifdef USE_GEMS3K_SERVER
#ifdef NO_ASYNC_SERVER
#include "gemsreaktoro/zmq_req_client.hpp"
#else
#include "gemsreaktoro/zmq_client.hpp"
#endif
#endif

TProfil* TProfil::pm;

const std::string _GEMS_version_stamp = std::string(" GEMS-GUI v.") + GEMSGUI_VERSION + " c." + GEMSGUI_VERSION_HASH;

// Heading lines of the Controls page (read-only fields, not part of the record)
static const char pa_heading_ipm[64] =    "IPM solver settings ---------------------------------------";
static const char pa_heading_optima[64] = "Optima solver settings ------------------------------------";

SPP_SETTING pa_ = {
    "Tolerances and controls: GEMSGUI v." GEMSGUI_VERSION " and GEMS3K v." GEMS3K_VERSION,
    {   // Typical default set (30.06.2025) new PSSC( logSI ) & uDD(), 0.7,   /* AG */   -0.0065,
        2,  /* PC */  2,     /* PD */   -4,   /* PRD */
        1,  /* PSM  */ 130,  /* DP */   1,   /* DW */
        0, /* DT */     30000,   /* PLLG */   1,  /* PE */  9999, /* IIM */
        1000., /* DG */   1e-13,  /* DHB */  1e-20,  /* DS */
        1e-5,  /* DK */  0.01,  /* DF */  0.01,  /* DFM */
        1e-5,  /* DFYw */  1e-5,  /* DFYaq */    1e-5,  /* DFYid */
        1e-5,  /* DFYr,*/  1e-5,  /* DFYh,*/   1e-5,  /* DFYc,*/
        1e-6, /* DFYs, */  1e-17,  /* DB */   0.7,   /* AG */
        -0.0065,   /* DGC */   1.0,   /* GAR */  1000., /* GAH */
        1e-3, /* GAS */   12.05,  /* DNS */   1e-13,  /* XwMin, */
        1e-13,  /* ScMin, */  1e-33, /* DcMin, */   1e-20, /* PhMin, */
        1e-5,  /* ICmin */   1e-10,  /* EPS */   1e-3,  /* IEPS */
        1e-10,  /* DKIN  */ nullptr,  /* tprn */
    },
    "CSC   ",   /* DCpct[6] */  "OjjbC+---",    /* DCpdc[9] */
    "+----M",   /* BCpc[6] */   "K3C   ",   /* REpct[6] */
    "SjjbC ",   /* REpdc[6]*/  "+-------", /* REpvc[8]  */
    "rjjbC-TCb-",   /* RPpdc[10]  */
    "++*-1-23----------4-------------", /* RPpvc[32] reserved free */
    "INNINN",   /* PHsol_t[6] */  "s-----",   /* PHpvc[6] */
    "++++-+-+++", /* MUpmv[10] */ "jjbC++0+", /* TPpdc[8] */
    "*-*-*-*-*-----------", /* TPpvc[20] */ "+-+-+-----", /* SYppc[10] */
    "*****-*-*-**-***-----------*", /* SYpvc[28]*/  "-+++----ME", /* UTppc[10] */
    "0*-------P--", /* PEpsc[12]  */  "------------", /* PEpvc[12] */
    { "GTDEMO task name   ", "Graphic screen # " } ,   /* GDcode[2][20] */
    "Plot ",                  /* GDpsc[7] */
    { "abscissa","ordinate"},    /* GDpcc[2][9] */
    "++++-+",   /* GDptc[6] */  "++++--",   /* GDpgw[6] */
    "*",   /* sdrefkey 32*/  "",   /* Reserv [50-32] */
    /* RTPARM */
    1, 15, 0,   /* P T V */  2,  /* Mode */
    { 1, 2, 3, 4, 5, },  /* ResShort */
    /* RTPARM */
    {0., 1., 0.}, /* Pi[3]*/  {0.0, 350., 25.}, /* Ti[3] */
    {0., 0., 0.},   /* Vi[3] */
    1., 25.,      /* Pr, Tr */
    /* MULTI */
    1.66e-24f, /* lowPosNum, */  1e-16f, /* logXw, */   1e-9f,  /* logYFk */
    { 0.064f, 3.72f, 0., 0., 0. },  // Defaults for aqueous auto-phase models
}; /* SPP_SETTING */

BASE_PARAM dfBase[5] =
        {     // Added to facilitate pre-setting in projects (SD,DK), revised 24.03.2020 DK
            { // Variant for aquatic systems with moderate non-ideality
                2,  /* PC */  2,     /* PD */   -4,   /* PRD */
                1,  /* PSM  */ 130,  /* DP */   1,   /* DW */
                0, /* DT */     30000,   /* PLLG */   1,  /* PE */
                9999,   /* IIM */
                1000., /* DG */   1e-13,  /* DHB */  1e-20,  /* DS */
                1e-5,  /* DK */  0.01,  /* DF */  0.01,  /* DFM */
                1e-5,  /* DFYw */  1e-5,  /* DFYaq */    1e-5,  /* DFYid */
                1e-5,  /* DFYr,*/  1e-5,  /* DFYh,*/   1e-5,  /* DFYc,*/
                1e-6, /* DFYs, */  1e-17,  /* DB */   1.,   /* AG */
                -0.98,   /* DGC */   1.0,   /* GAR */  1000., /* GAH */
                1e-3, /* GAS */   12.05,  /* DNS */   1e-13,  /* XwMin, */
                1e-13,  /* ScMin, */  1e-33, /* DcMin, */   1e-20, /* PhMin, */
                1e-5,  /* ICmin */   1e-10,  /* EPS */   1e-3,  /* IEPS */
                1e-10,  /* DKIN  */ nullptr,  /* tprn */
	    },
            { // Variant for aquatic systems with SCMs in sorption phases
                2,  /* PC */  2,     /* PD */   -4,   /* PRD */
                1,  /* PSM  */ 130,  /* DP */   1,   /* DW */
                0, /* DT */     -30000,   /* PLLG */   1,  /* PE */
                7000,   /* IIM */
                1000., /* DG */   1e-13,  /* DHB */  1e-20,  /* DS */
                3e-5,  /* DK */  0.01,  /* DF */  0.01,  /* DFM */
                1e-5,  /* DFYw */  1e-5,  /* DFYaq */    1e-5,  /* DFYid */
                1e-5,  /* DFYr,*/  1e-5,  /* DFYh,*/   1e-5,  /* DFYc,*/
                1e-6, /* DFYs, */  1e-17,  /* DB */   1.0,   /* AG */
                -0.01,   /* DGC */   1.0,   /* GAR */  1000., /* GAH */
                1e-3, /* GAS */   12.05,  /* DNS */   1e-13,  /* XwMin, */
                1e-13,  /* ScMin, */  1e-33, /* DcMin, */   1e-20, /* PhMin, */
                1e-5,  /* ICmin */   1e-10,  /* EPS */   5e-3,  /* IEPS */
                1e-10,  /* DKIN  */ nullptr,  /* tprn */
	    },
            { // Variant for fluid-rock systems with highly non-ideal phases
                2,  /* PC */  2,     /* PD */   -4,   /* PRD */
                1,  /* PSM  */ 130,  /* DP */   1,   /* DW */
                1, /* DT */     30000,   /* PLLG */   1,  /* PE */
                9990,   /* IIM */
                1000., /* DG */   1e-11,  /* DHB */  1e-20,  /* DS */
                6e-5,  /* DK */  0.01,  /* DF */  0.001,  /* DFM */
                1e-5,  /* DFYw */  1e-5,  /* DFYaq */    1e-5,  /* DFYid */
                1e-5,  /* DFYr,*/  1e-5,  /* DFYh,*/   1e-5,  /* DFYc,*/
                1e-6, /* DFYs, */  1e-17,  /* DB */   1.,   /* AG */
                -0.001,   /* DGC */   1.0,   /* GAR */  1000., /* GAH */
                1e-3, /* GAS */   12.05,  /* DNS */   1e-11,  /* XwMin, */
                1e-11,  /* ScMin, */  1e-30, /* DcMin, */   1e-20, /* PhMin, */
                1e-5,  /* ICmin */   1e-10,  /* EPS */   1e-3,  /* IEPS */
                1e-10,  /* DKIN  */ nullptr,  /* tprn */
            },
            { // Variant of strict settings for highest accuracy
                2,  /* PC */  2,     /* PD */   -5,   /* PRD */
                1,  /* PSM  */ 130,  /* DP */   1,   /* DW */
                0, /* DT */     30000,   /* PLLG */   1,  /* PE */
                7000,   /* IIM */
                1000., /* DG */   1e-14,  /* DHB */  1e-20,  /* DS */
                1e-6,  /* DK */  0.01,  /* DF */  0.001,  /* DFM */
                1e-5,  /* DFYw */  1e-5,  /* DFYaq */    1e-5,  /* DFYid */
                1e-5,  /* DFYr,*/  1e-5,  /* DFYh,*/   1e-5,  /* DFYc,*/
                1e-6, /* DFYs, */  1e-17,  /* DB */   1.,   /* AG */
                0.,   /* DGC */   1.0,   /* GAR */  1000., /* GAH */
                1e-3, /* GAS */   12.05,  /* DNS */   1e-13,  /* XwMin, */
                1e-13,  /* ScMin, */  1e-40, /* DcMin, */   1e-30, /* PhMin, */
                1e-7,  /* ICmin */   1e-10,  /* EPS */   1e-3,  /* IEPS */
                1e-10,  /* DKIN  */ nullptr,  /* tprn */
            },
            { // Variant of strict settings for reactive mass transport (with sorption)
                2,  /* PC */  2,     /* PD */   -4,   /* PRD */
                1,  /* PSM  */ 130,  /* DP */   1,   /* DW */
                0, /* DT */     30000,   /* PLLG */   1,  /* PE */
                7000,   /* IIM */
                1000., /* DG */   1e-14,  /* DHB */  1e-20,  /* DS */
                1e-5,  /* DK */  0.01,  /* DF */  0.01,  /* DFM */
                1e-5,  /* DFYw */  1e-5,  /* DFYaq */    1e-5,  /* DFYid */
                1e-5,  /* DFYr,*/  1e-5,  /* DFYh,*/   1e-5,  /* DFYc,*/
                1e-6, /* DFYs, */  1e-17,  /* DB */   1.,   /* AG */
                -0.001,   /* DGC */   1.0,   /* GAR */  1000., /* GAH */
                1e-3, /* GAS */   12.05,  /* DNS */   1e-13,  /* XwMin, */
                1e-13,  /* ScMin, */  1e-33, /* DcMin, */   1e-20, /* PhMin, */
                1e-5,  /* ICmin */   1e-10,  /* EPS */   1e-3,  /* IEPS */
                1e-10,  /* DKIN  */ nullptr,  /* tprn */
            }
	    // More preset parameter structures can be added here
            // (if so, increase number of elements from 5 above and in dfBase() call below)
	};

// Setup one of 5 default IPM numerical settings
void TProfil::ChangeSettings(int nSettings)
{
	if( nSettings == 0  )
     return; // No change in controls set
	else if( nSettings == 1 )
		    pa.p = pa_.p;
	     else
            pa.p = dfBase[std::min(nSettings-2,5)];
    // Added by DK 16.02.2012
    std::string pa_ver = _GEMS_version_stamp;
    pa_ver += _GEMIPM_version_stamp;
    memcpy(pa.ver, pa_ver.c_str(), TDBVERSION-1);
    pa.ver[TDBVERSION-1]='\0';
}

// Resets the BASE_PARAM fields stored in the tail of the project record to their defaults.
// Needed before reading a project saved by an older version, whose record lacks these objects.
void TProfil::resetExtendedSettings()
{
    pa.p.PSTALL = pa_.p.PSTALL;
    pa.p.IpmStallWindow = pa_.p.IpmStallWindow;
    pa.p.MbReproject = pa_.p.MbReproject;
    pa.p.DeterminacyWarn = pa_.p.DeterminacyWarn;
    pa.p.ColdRetryNudges = pa_.p.ColdRetryNudges;
    pa.p.FilloutBudget = pa_.p.FilloutBudget;
    pa.p.StabTPD = pa_.p.StabTPD;
    pa.p.IpmAugmentedKKT = pa_.p.IpmAugmentedKKT;
    pa.p.OptimaTol = pa_.p.OptimaTol;
    pa.p.LogBarrierTau = pa_.p.LogBarrierTau;
    pa.p.PhaseHessianFloor = pa_.p.PhaseHessianFloor;
    pa.p.OptimaStallWindow = pa_.p.OptimaStallWindow;
    pa.p.OptimaMaxSeconds = pa_.p.OptimaMaxSeconds;
    pa.p.OptimaFDHessian = pa_.p.OptimaFDHessian;
    pa.p.OptimaMoleFracHessian = pa_.p.OptimaMoleFracHessian;
    pa.p.OptimaFDHessianDelay = pa_.p.OptimaFDHessianDelay;
    pa.p.OptimaDcFloor = pa_.p.OptimaDcFloor;
    pa.p.OptimaEarlyStabilityAt = pa_.p.OptimaEarlyStabilityAt;
    pa.p.OptimaDimReduce = pa_.p.OptimaDimReduce;
    pa.p.OptimaDimReduceTol = pa_.p.OptimaDimReduceTol;
    pa.p.OptimaPreSolveFirstIters = pa_.p.OptimaPreSolveFirstIters;
    pa.p.OptimaZeroAbsent = pa_.p.OptimaZeroAbsent;
    pa.p.OptimaLineSearch = pa_.p.OptimaLineSearch;
    pa.p.OptimaLSStallEscape = pa_.p.OptimaLSStallEscape;
    pa.p.OptimaLSRejectWorse = pa_.p.OptimaLSRejectWorse;
    pa.p.OptimaTpdAccept = pa_.p.OptimaTpdAccept;
    pa.p.OptimaCgSeed = pa_.p.OptimaCgSeed;
    pa.p.OptimaColdRetry = pa_.p.OptimaColdRetry;
    pa.p.OptimaFinish = pa_.p.OptimaFinish;
    pa.p.OptimaAcceptRepair = pa_.p.OptimaAcceptRepair;
    pa.p.MbClassRule = pa_.p.MbClassRule;
    pa.p.FailRescue = pa_.p.FailRescue;
    pa.p.OptimaRedoxFree = pa_.p.OptimaRedoxFree;
    pa.p.OptimaDryOut = pa_.p.OptimaDryOut;
    pa.p.OptimaCallCap = pa_.p.OptimaCallCap;
}

void BASE_PARAM::write(GemDataStream& oss)
{
    oss.writeArray( &PC, 10 );
    oss.writeArray( &DG, 28 );
    oss.writeArray( "0000", 4 );
}

void BASE_PARAM::read(GemDataStream& iss)
{
    char tmp[4];
    iss.readArray( &PC, 10 );
    iss.readArray( &DG, 28 );
    iss.readArray( tmp, 4 );
}


void SPP_SETTING::write(GemDataStream& oss)
{
    oss.writeArray( ver, TDBVERSION );
    p.write( oss );
    oss.writeArray( DCpct, 352 );
    oss.writeArray( &NP, 9 );
    oss.writeArray( Pi, 19 );
}

void SPP_SETTING::read(GemDataStream& iss)
{
    iss.readArray( ver, TDBVERSION );
    p.read( iss );
    iss.readArray( DCpct, 352 );
    iss.readArray( &NP, 9 );
    iss.readArray( Pi, 19 );
}


TProfil::TProfil( uint nrt ):
        TCModule( nrt )
{
    aFldKeysHelp.push_back("Name of the modeling project");
    aFldKeysHelp.push_back("Comment to the project definition");
    //startKeyEdit = 0;
    start_title = " Numerical and Configuration Settings ";
    pa= pa_;
    pa.p.tprn=nullptr;

    rmults = nullptr;
    mtparm =  nullptr;
    syst = nullptr;
    multi_internal =nullptr;

    //mup = 0;
    //tpp = 0;
    //syp = 0;
    //pmp = 0;

    SFold = nullptr;
    SMold = nullptr;
    SAold = nullptr;
    SBold = nullptr;
    Llold =nullptr;

    userCancel1 = false;
    stepWise1 = false;
    calcFinished = false;
    comp_change_all = false;
    internalBufer = nullptr;
}

// init submodules to calc module
void TProfil::InitSubModules()
{
    if( !rmults /*aMod.GetCount() < MD_RMULTS*/ )
    {
        rmults = new TRMults( MD_RMULTS );
        aMod.push_back( std::shared_ptr<TRMults>(rmults) );
        TRMults::sm = rmults;
        rmults->ods_link();
        //mup = rmults->GetMU();

        aMod.push_back( std::shared_ptr<TMTparm>( mtparm = new TMTparm( MD_MTPARM ) ));
        TMTparm::sm = mtparm;
        mtparm->ods_link();
        //tpp = mtparm->GetTP();
        aMod.push_back( std::shared_ptr<TSyst>( syst = new TSyst( MD_SYSTEM )) );
        TSyst::sm = syst;
        syst->ods_link();
        //syp = syst->GetSY();
        aMod.push_back( std::shared_ptr<TMulti>(multi_internal = new TMulti( MD_MULTI )) );
        TMulti::sm = multi_internal;
        //multi_internal->setPa(this);
        //pmulti = multi;
        multi_internal->ods_link();
        //pmp = multi->GetPM();
        aMod.push_back( std::shared_ptr<TEQCalc>( new TEQCalc( MD_EQCALC )) );
        aMod.push_back( std::shared_ptr<TEQDemo>(new TEQDemo( MD_EQDEMO )) );

    }
}

// link values to objects
void TProfil::ods_link( int )
{
    /* rmults->ods_link();
     mtparm->ods_link();
     syst->ods_link();
     multi->ods_link();*/
    aObj[o_paver]->SetPtr( pa.ver );
    aObj[o_paxhd1]->SetPtr( const_cast<char*>(pa_heading_ipm) );
    aObj[o_paxhd2]->SetPtr( const_cast<char*>(pa_heading_optima) );
    aObj[o_papc]->SetPtr( &(pa.p.PC) );
    aObj[o_paprd]->SetPtr(&(pa.p.PD) );
    aObj[o_padpwt]->SetPtr( &(pa.p.DP));
    aObj[o_papllg]->SetPtr( &(pa.p.PLLG));
    aObj[o_pape]->SetPtr( &(pa.p.PE) );
    aObj[o_paiim]->SetPtr( &(pa.p.IIM) );
    aObj[o_padg]->SetPtr( &(pa.p.DG) );
    if( pa.p.DG <1 )                     // SD 07/05/2010
        pa.p.DG = 1e4;
    aObj[o_padhb]->SetPtr( &(pa.p.DHB) );
    aObj[o_pads]->SetPtr( &(pa.p.DS) );
    aObj[o_padk]->SetPtr( &(pa.p.DK) );
    aObj[ o_padf]->SetPtr(   &(pa.p.DF));
    aObj[ o_padfy]->SetPtr(  &(pa.p.DFYw));
    aObj[ o_padb]->SetPtr(   &(pa.p.DB));
    aObj[ o_paag]->SetPtr(   &(pa.p.AG));
    aObj[ o_padgc]->SetPtr(  &(pa.p.DGC));
    aObj[ o_pagan]->SetPtr(  &(pa.p.GAR));
    aObj[ o_padns]->SetPtr(  &(pa.p.DNS));
    aObj[ o_paxmin]->SetPtr( &(pa.p.XwMin));
    aObj[ o_paeps]->SetPtr(  &(pa.p.EPS));
    aObj[ o_padkin]->SetPtr( &(pa.p.DKIN));
    aObj[ o_patprn]->SetPtr( pa.p.tprn );
    aObj[ o_padcpct]->SetPtr(pa.DCpct );
    aObj[ o_padcpdc]->SetPtr(pa.DCpdc );
    aObj[ o_parepct]->SetPtr(pa.REpct );
    aObj[ o_parepdc]->SetPtr(pa.REpdc );
    aObj[ o_parepvc]->SetPtr(pa.REpvc );
    aObj[ o_parppdc]->SetPtr(pa.RPpdc );
    aObj[ o_parppvc]->SetPtr(pa.RPpvc );
    aObj[ o_paphsol_t]->SetPtr(pa.PHsol_t );
    aObj[ o_paphpvc]->SetPtr(pa.PHpvc );
    aObj[ o_pamupmv]->SetPtr(pa.MUpmv );
    aObj[ o_patppdc]->SetPtr(pa.TPpdc );
    aObj[ o_patppvc]->SetPtr(pa.TPpvc );
    aObj[ o_pasyppc]->SetPtr(pa.SYppc );
    aObj[ o_pasypvc]->SetPtr(pa.SYpvc );
    aObj[ o_pautppc]->SetPtr(pa.UTppc );
    aObj[ o_papepsc]->SetPtr(pa.PEpsc );
    aObj[ o_papepvc]->SetPtr(pa.PEpvc );
    aObj[ o_pagdcode]->SetPtr(&pa.GDcode[0][0] );
    aObj[ o_pagdpsc]->SetPtr(pa.GDpsc );
    aObj[ o_pagdpcc]->SetPtr(&pa.GDpcc[0][0] );
    aObj[ o_pagdptc]->SetPtr(pa.GDptc );
    aObj[ o_parpnptv]->SetPtr(&pa.NP );
    aObj[ o_parpmode]->SetPtr(&pa.Mode );
    aObj[ o_partpi]->SetPtr(   pa.Pi );
    aObj[ o_partti]->SetPtr(   pa.Ti );
    aObj[ o_partvi]->SetPtr(   pa.Vi );
    aObj[ o_padrprtr]->SetPtr(&pa.DRpst );
    aObj[ o_papmnum]->SetPtr( &pa.lowPosNum );

    aObj[ o_pabcpc]->SetPtr(  pa.BCpc );
    aObj[ o_pagdpgw]->SetPtr(  pa.GDpgw );
    aObj[ o_pasdref]->SetPtr(  pa.SDrefKey );

    aObj[ o_paxstall]->SetPtr( &(pa.p.PSTALL));
    aObj[ o_paxipmsw]->SetPtr( &(pa.p.IpmStallWindow));
    aObj[ o_paxmbrep]->SetPtr( &(pa.p.MbReproject));
    aObj[ o_paxdetw]->SetPtr( &(pa.p.DeterminacyWarn));
    aObj[ o_paxcrn]->SetPtr( &(pa.p.ColdRetryNudges));
    aObj[ o_paxfob]->SetPtr( &(pa.p.FilloutBudget));
    aObj[ o_paxstpd]->SetPtr( &(pa.p.StabTPD));
    aObj[ o_paxkkt]->SetPtr( &(pa.p.IpmAugmentedKKT));
    aObj[ o_paxotol]->SetPtr( &(pa.p.OptimaTol));
    aObj[ o_paxlbt]->SetPtr( &(pa.p.LogBarrierTau));
    aObj[ o_paxphf]->SetPtr( &(pa.p.PhaseHessianFloor));
    aObj[ o_paxosw]->SetPtr( &(pa.p.OptimaStallWindow));
    aObj[ o_paxoms]->SetPtr( &(pa.p.OptimaMaxSeconds));
    aObj[ o_paxofh]->SetPtr( &(pa.p.OptimaFDHessian));
    aObj[ o_paxomh]->SetPtr( &(pa.p.OptimaMoleFracHessian));
    aObj[ o_paxofd]->SetPtr( &(pa.p.OptimaFDHessianDelay));
    aObj[ o_paxodf]->SetPtr( &(pa.p.OptimaDcFloor));
    aObj[ o_paxoes]->SetPtr( &(pa.p.OptimaEarlyStabilityAt));
    aObj[ o_paxodr]->SetPtr( &(pa.p.OptimaDimReduce));
    aObj[ o_paxodt]->SetPtr( &(pa.p.OptimaDimReduceTol));
    aObj[ o_paxops]->SetPtr( &(pa.p.OptimaPreSolveFirstIters));
    aObj[ o_paxoza]->SetPtr( &(pa.p.OptimaZeroAbsent));
    aObj[ o_paxols]->SetPtr( &(pa.p.OptimaLineSearch));
    aObj[ o_paxole]->SetPtr( &(pa.p.OptimaLSStallEscape));
    aObj[ o_paxorw]->SetPtr( &(pa.p.OptimaLSRejectWorse));
    aObj[ o_paxota]->SetPtr( &(pa.p.OptimaTpdAccept));
    aObj[ o_paxocs]->SetPtr( &(pa.p.OptimaCgSeed));
    aObj[ o_paxocr]->SetPtr( &(pa.p.OptimaColdRetry));
    aObj[ o_paxofn]->SetPtr( &(pa.p.OptimaFinish));
    aObj[ o_paxoar]->SetPtr( &(pa.p.OptimaAcceptRepair));
    aObj[ o_paxmcr]->SetPtr( &(pa.p.MbClassRule));
    aObj[ o_paxfrs]->SetPtr( &(pa.p.FailRescue));
    aObj[ o_paxorf]->SetPtr( &(pa.p.OptimaRedoxFree));
    aObj[ o_paxodo]->SetPtr( &(pa.p.OptimaDryOut));
    aObj[ o_paxocc]->SetPtr( &(pa.p.OptimaCallCap));
    aObj[ o_spppar]->SetPtr(  static_cast<void *>(&pa) );
    aObj[ o_spppar]->SetM( sizeof( SPP_SETTING ) );
aObj[ o_sptext]->SetPtr(  internalBufer );
    //   aObj[ o_sppconst]->SetPtr( sc_ );
    //   ob[o_sppconst].dim_M = sizeof( SPP_CONST );
    //   aObj[ o_sppdatac]->SetPtr( &nQ );
    //   ob[o_sppdatac].dim_M = sizeof( DATACOUNT );
}

// set dynamic Objects ptr to values

void TProfil::dyn_set(int )
{
    pa.p.tprn= static_cast<char *>(aObj[o_patprn]->GetPtr());
    internalBufer = static_cast<char *>(aObj[ o_sptext]->GetPtr());
    charge_mismatch_quest_reply = VF_UNDEF;
    if( rmults ) rmults->dyn_set();
    if( mtparm ) mtparm->dyn_set();
    if( syst ) syst->dyn_set();
    if( multi_internal ) multi_internal->dyn_set();
}

// free dynamic memory in objects and values
void TProfil::dyn_kill(int )
{
    pa.p.tprn = static_cast<char *>(aObj[o_patprn]->Free());
internalBufer = static_cast<char *>(aObj[ o_sptext]->Free());
    if( rmults ) rmults->dyn_kill();
    if( mtparm ) mtparm->dyn_kill();
    if( syst ) syst->dyn_kill();
    if( multi_internal ) multi_internal->dyn_kill();
}

// realloc dynamic memory
void TProfil::dyn_new(int )
{
    /*  rmults->dyn_new();
      mtparm->dyn_new();
      syst->dyn_new();
      multi->dyn_new();
    */
    // if( pa.p.tprn == 0 )
    //   pa.p.tprn = (char *)aObj[o_patprn].Alloc( 1, 256, S_ );
}


//set default information

void TProfil::set_def( int )
{
    pa = pa_;
    internalBufer = nullptr;
    charge_mismatch_quest_reply = VF_UNDEF;
    if( rmults ) rmults->set_def();
    if( mtparm ) mtparm->set_def();
    if( syst ) syst->set_def();
    if( multi_internal ) multi_internal->set_def();
}

/* opens window with 'Remake record' parameters
*/
void TProfil::MakeQuery()
{
    const char * p_key;
    char flgs[38];
    int tasktype = 0;

    p_key  = db->PackKey();
    // flgs 0-9 MSpmv, 10-29 TPptv, 30-33  TPun, 34-37 TPsv
    memcpy( flgs, &rmults->GetMU()->PmvSA, 10);
    memcpy( flgs+10, &mtparm->GetTP()->PtvG, 20);
    memcpy( flgs+30, &mtparm->GetTP()->PunE, 8);

    if( !vfProjectSet( pVisor->window()/*window()*/, p_key, flgs, tasktype ))
         Error( p_key, "Project record configuration cancelled by the user!" );

    memcpy( &rmults->GetMU()->PmvSA, flgs,    10);
    memcpy( &mtparm->GetTP()->PtvG,  flgs+10, 20);
    memcpy( &mtparm->GetTP()->PunE,  flgs+30, 8);
    // SD 07/08/2009 set BASE_PARAM from default
    ChangeSettings(tasktype);
}

// Help on Modelling Project module ( ? button )
const char* TProfil::GetHtml()
{
   if(pVisor->ProfileMode )
     return GEMS_NUMSET_HTML;
   else
    return GM_PROJECT_HTML;
}

void TProfil::makeGEM2MTFiles(QWidget* par )
{
    MULTI *pmp = multi_internal->GetPM();

    try
	 {
      // set default data and realloc arrays
      char flags[6];
      double Tai[4], Pai[4];
      long int nTp_, nPp_;
      double *arT, *arP;

      Tai[0] = pmp->TCc; Tai[1] = pmp->TCc;
      Pai[0] = pmp->Pc; Pai[1] = pmp->Pc;

      if( !vfLookupDialogSet(par, flags, Tai, Pai ) )// ask constants
        	return;

      //internal objects for lookup arrays data
      arT = static_cast<double *>(aObj[ o_w_tval]->GetPtr());
      nTp_ = aObj[ o_w_tval]->GetN();
      arP = static_cast<double *>(aObj[ o_w_pval]->GetPtr());
      nPp_ = aObj[ o_w_pval]->GetN();

      auto na = TNodeArrayGUI::create(1, multi_internal);
      na->setExportNode();

      // realloc and setup data for dataCH and DataBr structures
      na->MakeNodeStructuresOne( par, ( flags[0] == S_OFF ),( flags[4] == S_ON ),
          arT, arP, nTp_, nPp_, Tai[3], Pai[3]  );

       // make  all files
      if(flags[5] == S_OFF )
           na->pNodT0()[0]->NodeStatusFMT = No_nodearray;

      na->PutGEM2MTFiles( par, 1,  flags[1], ( flags[2] == S_ON ),
                   ( flags[3] == S_ON ), false, true );// addMui, to txt
      //Test to compare standalone  na->GEM_print_ipm( "GemsCalcPhase.txt" );
      //outMultiTxt( "IPM_GEMS.txt"  );
    }
    catch( TError& xcpt )
    {
      throw;
    }

}

// Checks if the modelling project name is the same as read from GEMS3K I/O files
bool TProfil::CompareProjectName( const char* SysKey )
{
    auto len = rt[RT_PARAM]->FldLen(0);
    const char* proj_key = db->UnpackKey();
    if( memcmp( SysKey, proj_key, len ) )
        return true;
    else
      return false;
}

// Change the key read from GEMS3K I/O files
// Copy T and P from DATABR
void TProfil::ChangeTPinKey( double T, double P )
{
    MULTI* pmp = multi_internal->GetPM();
    char bT[40];
    char bP[40];

    Gcvt( T, 6, bT );
    Gcvt( P, 6, bP );
    rt[RT_SYSEQ]->SetKey( pmp->stkey );
    rt[RT_SYSEQ]->MakeKey( RT_SYSEQ, pmp->stkey, RT_SYSEQ, 0, RT_SYSEQ, 1,
                           RT_SYSEQ, 2, RT_SYSEQ, 3, RT_SYSEQ, 4,
                           K_IMM, bP, K_IMM, bT, RT_SYSEQ, 7, K_END);

    std::string str = pmp->stkey;
    std::string capName = "Change the key read from GEMS3K I/O files";
AGAIN:
    str = TSysEq::pm->GetKeyofRecord( str.c_str(),capName.c_str(), KEY_NEW );
    if( str.empty() )
         Error( GetName(), "Record creation rejected!");
    if( rt[RT_SYSEQ]->FindCurrent( str.c_str() ) >= 0 )
    {
        capName = "This record already exists! Please, enter another name.";
        goto AGAIN;
    }
    rt[RT_SYSEQ]->SetKey( str.c_str() );
}




//  This function sets the system/SysEq switches for components and phases
//       according to mui, muj, mup index lists in MULTI that were read in.
void TProfil::SetSysSwitchesFromMulti( )
{
     MULTI* pmp = multi_internal->GetPM();
     RMULTS* mup = rmults->GetMU();
     SYSTEM* syp = syst->GetSY();
     long i, ii, j, jj, k, kk;


     // Set default informations to arrays
     syst->setDefData();

     // ICs
     for( ii=0; ii < mup->N; ii++ )
        syp->Icl[ii] = S_OFF;
     for( i=0; i < pmp->N; i++ )
     {
        ii = pmp->mui[i];
        syp->Icl[ii] = S_ON;
     }
     syp->N = pmp->N;
     // DCs
     for( jj=0; jj < mup->L; jj++ )
     {
         syp->Dcl[jj] = S_OFF;
     //    syp->DLL[jj] = 0.;
     //    syp->DUL[jj] = 1e6;
     }
     for( j=0; j < pmp->L; j++ )
     {
        jj = pmp->muj[j];
        syp->Dcl[jj] = S_ON;
     //   syp->DLL[jj] = pmp->DLL[j];
     //   syp->DUL[jj] = pmp->DUL[j];
     }
     syp->L = pmp->L;
     syp->Ls = pmp->Ls;
     // Phases
     for( kk=0; kk < mup->Fi; kk++ )
        syp->Pcl[kk] = S_OFF;
     for( k=0; k < pmp->FI; k++ )
     {
        kk = pmp->muk[k];
        syp->Pcl[kk] = S_ON;
        syp->Ll[kk] = pmp->L1[k];
     }
     syp->Fi = pmp->FI;
     syp->Fis = pmp->FIs;
}

// Reading structure MULTI (GEM IPM work structure)
void TProfil::CmReadMulti(const char* path)
{
    MULTI* pmp = multi_internal->GetPM();
    SYSTEM* syp = syst->GetSY();
    TNodeGUI na(multi_internal);

    if( na.GEM_init( path ) )
    {
      Error( path, "GEMS3K Init() error: \n" + na.ipmLogError() );
    }
    multi_internal->dyn_set();

    // Here to compare the modelling project name; error when from a different project.
    if( CompareProjectName( pmp->stkey ) )
    {
        Error( pmp->stkey, "E15IPM: Wrong project name by reading GEMS3K I/O files ");
    }

     // Unpacking the actual contents of DBR file including speciation
    na.unpackDataBr( true );
    for( int j=0; j < pmp->L; j++ )
        pmp->X[j] = pmp->Y[j];
    pmp->TC = pmp->TCc;
    pmp->T =  pmp->Tc;
    pmp->P =  pmp->Pc;
#ifndef USE_GEMS3K_SERVER
    pmp->VX_ = pmp->VXc; // from cm3 to m3
#endif
    // Set T and P  for key from DataBr
    ChangeTPinKey( pmp->TC, pmp->P );

    pmp->pESU = 2;  // SysEq unpack flag set

    // Test sizes and flags
    // my be diffirent compared to internal system (do not test):
    // syp->PAalp, syp->PSigm, pmp->FIat, pmp->E
    int iRet = BAL_compare();
    switch(iRet)
    {
      case 0: // Some lists have changed in the loaded GEM task relative to base system
           pmp->pIPN = 0;
           if( pmp->pTPD > 1)
               pmp->pTPD = 1; // reload Go, Vol
           pmp->pKMM = 0;
           // sets the system/SysEq switches for
           // components and phases according to mui, muj, mup that were read in.
           SetSysSwitchesFromMulti( );
//           vfMessage(par, "Changes in the GEMS3K system setup relative to GEMS",
//                    "Some vectors were allocated or some dimensions changed.\n"
//                    "The equilibrium cannot be recalculated");
          break;
       case 1: // Bulk composition has changed, but lists are the same
          //vfMessage(par, "Changes in the GEMS3K system setup relative to GEMS",
          //          "Some constraints are not the same.\n"
          //          "Different results may be obtained upon recalculation");
       default:
          break;
    }

    // for loading GEX to System
    CheckMtparam();
    multi_internal->TMultiBase::DC_LoadThermodynamicData(&na);

    // Unpack the pmp->B vector (b) into syp->BI and syp->BI (BI_ vector).

    for( long i=0; i < pmp->N; i++ ) {
#ifndef  USE_GEMS3K_SERVER
        syp->BI[pmp->mui[i]] =
#endif
        syp->B[pmp->mui[i]] = pmp->B[i];
    }

    for( long jj, j=0; j < pmp->L; j++ )
    {
       jj = pmp->muj[j];
       syp->DLL[jj] = pmp->DLL[j];
       syp->DUL[jj] = pmp->DUL[j];
       syp->GEX[jj] = na.DC_G0(j, pmp->P*bar_to_Pa, pmp->Tc, false) - mtparm->GetTP()->G[jj];
    }

// !!! We cannot restore to System  for adsorbtion must be done

    TSysEq::pm->CellChanged();

    // calculate mass of the system
     pmp->MBX = 0.0;
     for(int i=0; i<pmp->N; i++ )
           pmp->MBX += pmp->B[i] * pmp->Awt[i];
     pmp->MBX /= 1000.;

    // Restoring the rest of MULTI contents from primal and dual solution
    pmp->pIPN =0;
    multi_internal->Alloc_internal();
    multi_internal->EqstatExpand( /*pmp->stkey,*/ false );
    // outMultiTxt( "IPM_EqstatExpand.txt"  );
    // multi->Free_internal();
    // na.unpackDataBr( true );
    // We can get different results in GEMS than in GEMS3K
    // because of slightly different values into G, V ...
    // (interpolation of thermodynamic data or precision )
}



//Delete record with key
void
TProfil::DeleteRecord( const char *key, bool /*errifNo*/ )
{
    TCStringArray aList;
    TCIntArray anR;
    char pkey[81];
    uint i;

    int  Rnum = rt[nRT]->Find( key );
    ErrorIf( Rnum<0, GetName(), "Record to delete not found!");
    rt[RT_PARAM]->Get( Rnum ); // read record
    dyn_set();
    SetFN();                  // reopen files of data base
    rt[nRT]->SetKey( key);

    // Delete all records connected to project
    aList.clear();    //SYSEQ
    anR.clear();
    rt[RT_SYSEQ]->MakeKey( RT_PARAM, pkey, RT_PARAM, 0,
                           K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_END);
    rt[RT_SYSEQ]->GetKeyList( pkey, aList, anR );
    for( i=0; i< aList.size(); i++)
        TSysEq::pm->DeleteRecord(aList[i].c_str());

    aList.clear();    //PROCES
    anR.clear();
    rt[RT_PROCES]->MakeKey( RT_PARAM, pkey, RT_PARAM, 0,
                            K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_END);
    rt[RT_PROCES]->GetKeyList( pkey, aList, anR );
    for( i=0; i< aList.size(); i++)
        TProcess::pm->DeleteRecord(aList[i].c_str());

    aList.clear();    //UNSPACE
    anR.clear();
    rt[RT_UNSPACE]->MakeKey( RT_PARAM, pkey, RT_PARAM, 0,
      K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_END);
    rt[RT_UNSPACE]->GetKeyList( pkey, aList, anR );
    for( i=0; i< aList.size(); i++)
        TUnSpace::pm->DeleteRecord(aList[i].c_str());

    aList.clear();    //GTDEMO
    anR.clear();
    rt[RT_GTDEMO]->MakeKey( RT_PARAM, pkey, RT_PARAM, 0,
                            K_ANY, K_ANY, K_ANY, K_ANY, K_END);
    rt[RT_GTDEMO]->GetKeyList( pkey, aList, anR );
    for( i=0; i< aList.size(); i++)
        TGtDemo::pm->DeleteRecord(aList[i].c_str());

    aList.clear();    //DUALTH
    anR.clear();
    rt[RT_DUALTH]->MakeKey( RT_PARAM, pkey, RT_PARAM, 0,
                            K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_ANY, K_END);
    rt[RT_DUALTH]->GetKeyList( pkey, aList, anR );
    for( i=0; i< aList.size(); i++)
        TDualTh::pm->DeleteRecord(aList[i].c_str());

    rt[nRT]->Del( Rnum );
}


// Load Thermodynamic data from Database
void TProfil::CheckMtparam()
{

  if( fabs( mtparm->GetTP()->T - multi_internal->GetPM()->TCc ) > 1.e-10 ||
         fabs( mtparm->GetTP()->P - multi_internal->GetPM()->Pc ) > 1.e-9 )
   { // load new MTPARM on T or P
      mtparm->LoadMtparm( multi_internal->GetPM()->TCc, multi_internal->GetPM()->Pc );
      multi_internal->GetPM()->pTPD = 0;
   }
 }

// Setup of flags for MULTY remake
// pNP,  //Mode of FIA selection: 0-LPP AIA, 1-old eqstate, -1-user's choice
// pESU, // Unpack old eqstate from EQSTAT record?  0-no 1-yes
// pIPN, // State of IPN-arrays:  0-create; 1-available; -1 remake
// pBAL, // State of reloading CSD:  1- BAL only; 0-whole CSD
// pTPD, // State of reloading thermod data: 0- all    2 - no
void TProfil::PMtest( const char *key )
{
    //double V, T, P;
    TSysEq* STat = dynamic_cast<TSysEq *>(aMod[RT_SYSEQ].get());
    TProcess* Proc = dynamic_cast<TProcess *>(aMod[RT_PROCES].get());
    MULTI *pmp = multi_internal->GetPM();

    // test for available old solution
    if( STat->ifCalcFlag() )
    {
       if( !pmp->pESU )      // if pESU == 2 (task loaded before), left value
          pmp->pESU = 1;
    }
    else pmp->pESU = 0;

   // no old solution => must be LPP AIA
   if( pmp->pESU == 0 )
        pmp->pNP = 0;

   // test changes in the modified system relative to MULTI
   pmp->pBAL =  BAL_compare();
   gui_logger->debug("BAL_compare: pmp->pBAL {} pmp->pTPD {}", pmp->pBAL, pmp->pTPD);

   if( !pmp->pBAL ) // if some vectors were allocated or some dimensions changed
   {
       pmp->pIPN = 0;
       if( pmp->pTPD > 1)
          pmp->pTPD = 1; // reload Go, Vol
   }

    // special setup from process simulator
    if( Proc->pep->Istat == P_EXECUTE ||
        Proc->pep->Istat == P_MT_EXECUTE )
    {
        if(Proc->pep->PvR1 == S_OFF )
            pmp->pNP = 0;
        else
            pmp->pNP = 1;
    }

    multi_internal->MultiKeyInit( key );

}


// -------------------------------------------------------------------
// Compare changes in the modified system relative to MULTI
// if some vectors were allocated or some dimensions changed - return 0;
// else if bulk composition or some constraints has changed - return 1;
// else if nothing has changed - return 2.
short TProfil::BAL_compare()
{
    long i,j,k, jj, jb, je=0;
    double Go, Gg, Ge, pGo;
    SYSTEM *syp = TSyst::sm->GetSY();
    MULTI *pmp = multi_internal->GetPM();

    // Changes in thermdynamic (Saved DComp, ReacDC or Phase records )
    if( pmp->pTPD == -1 )  // 16/11/2011 SD
      return 0;

// Test A - sizes and selectors
    if( pmp->N != syp->N || pmp->L != syp->L || pmp->Ls != syp->Ls || pmp->LO != syp->Lw
            || pmp->PG != syp->Lg // || pmp->E != syp->PE  // check of E added 08.02.2010 DK
            || pmp->PSOL != syp->Lhc || pmp->Lads != syp->Lsor
            || pmp->FI != syp->Fi || pmp->FIs != syp->Fis )
        return 0;
    if( !((syp->PE != S_OFF && pmp->E == 1) || (syp->PE == S_OFF && pmp->E == 0))  )
        return 0;

    if(( syp->DLLim == S_ON || syp->DULim == S_ON ) && pmp->PLIM != 1 )
        return 0;
    else if( syp->DLLim == S_OFF && syp->DULim == S_OFF && pmp->PLIM == 1 )
        return 0;
// test selectors
    for( i=0; i<pmp->N; i++ )
        if( syp->Icl[pmp->mui[i]] == S_OFF )
            return 0;
    for( j=0; j<pmp->L; j++ )
        if( syp->Dcl[pmp->muj[j]] == S_OFF )
            return 0;
    for( k=0; k<pmp->FI; k++ )
        if( syp->Pcl[pmp->muk[k]] == S_OFF )
            return 0;


  // test other settings for G0
  for( k = 0; k < pmp->FI; k++ )
  {
	jb = je;
	je += pmp->L1[k];
    for( j=jb; j<je; j++ )
    {
       Gg = Ge = 0.0;    //   This part had to be changed after integrating Ge into pmp->G0
       jj = pmp->muj[j]; //        DK    07.03.2008,  16.05.2008

       Go = mtparm->GetTP()->G[jj]; //  G0(T,P) value taken from MTPARM
       if( syp->Guns )  // This is used mainly in UnSpace calculations
           Gg = syp->Guns[jj];    // User-set increment to G0 from project system
       if( syp->GEX && syp->PGEX != S_OFF )   // User-set increment to G0 from project system
           Ge = syp->GEX[jj];     //now Ge is integrated into pmp->G0 (since 07.03.2008) DK
       pGo = multi_internal->ConvertGj_toUniformStandardState( Go+Gg+Ge, j, k );

       if( fabs( pGo - pmp->G0[j] )* pmp->RT >= 0.001 )
       {
           break;   // GEX or Guns has changed for this DC in the system definition
       }
    }  // j
    if( j < je )
       return 0;
  }  // k


  // lists of components didn't change
  // test B - recipes and constraints
      for( i=0; i<pmp->N; i++ )
          if( noZero( pmp->B[i] ) )
          {
//              if( fabs( syp->B[pmp->mui[i]] /pmp->B[i] - 1. ) > pa.p.DHB * 0.001   )
              if( fabs( log10( syp->B[pmp->mui[i]] /pmp->B[i] )) > 1e-13  ) // log10( pa.p.DB ) ) )
                  return 1;
          }
          else
          {
              if( noZero( syp->B[pmp->mui[i]] )   )
                    return 1;
          }

    je = 0;
      // test other settings for DCs
    for( k = 0; k < pmp->FI; k++ )
    {
          jb = je;
          je += pmp->L1[k];
      for( j=jb; j<je; j++ )
      {
         jj = pmp->muj[j]; //        DK    07.03.2008,  16.05.2008

          if(( syp->DLLim != S_OFF ) && pmp->PLIM == 1 )
  //            if( fabs( syp->DLL[jj] - pmp->DLL[j] ) >= 1e-19 )
                if( !approximatelyEqual( syp->DLL[jj], pmp->DLL[j])  )   //SD 22/01/2009
                  break;
          if(( syp->DULim != S_OFF ) && pmp->PLIM == 1 )
  //            if( fabs( syp->DUL[jj] - pmp->DUL[j] ) >= 1e-19 )
                if( !approximatelyEqual( syp->DUL[jj], pmp->DUL[j] ) )   //SD 22/01/2009
                break;
          if( syp->DULim != S_OFF || syp->DLLim != S_OFF )
          {  if( pmp->RLC[j] != syp->RLC[jj] )
               break;
             if( pmp->RSC[j] != syp->RSC[jj] )
               break;
          }
      }  // j
      if( j < je )
         return 1;
    }  // k

    for( k=0; k<pmp->FI; k++ )
    {
      auto kk = pmp->muk[k];
      if( syp->PSigm != S_OFF )
      {
         if( fabs( pmp->Sigw[k] - syp->Sigm[kk][0]) > 1e-19 )
           break;
         if( fabs( pmp->Sigg[k] - syp->Sigm[kk][1]) > 1e-19 )
           break;
      }
      if( syp->PAalp != S_OFF )
        if( fabs( pmp->Aalp[k] - syp->Aalp[kk]) > 1e-19 )
         break;
    }

    if( k < pmp->FI )
        return 1;

    // bulk chem. compos. and constraints unchanged
    return 2;
}

//// Clear syp->XeD, syp->XeA, syp->Phm and copy pmp->B to syp->BI
void TProfil::Clear_XeA_XeD_Phm_BIun()
{
    int i;
    RMULTS* mup = rmults->GetMU();
    SYSTEM *syp = syst->GetSY();

     // put Bb data to system  ( into B and bi_)
    if( syp->PbIC != S_OFF )
     for( i=0; i<mup->N; i++)
      {
        syp->BI[i] =  0; //pmp->B[];
        syp->BIun[i] =  QUAN_MOL;
      }
// set zeros to xd_ xp_ and xa_
    if( syp->PbDC != S_OFF )
        for( i=0; i < mup->L; i++)
              syp->XeD[i] = 0.;
    if( syp->PbAC != S_OFF && rmults->GetMU()->La > 0 )
        for( i=0; i < mup->La; i++)
               syp->XeA[i] = 0.;
    if( syp->PbPH != S_OFF )
        for( i=0; i < mup->Fi; i++)
            syp->Phm[i] = 0.;
}


void TProfil::LoadFromMtparm(DATACH *CSD , bool no_interpolat)
{
    if( no_interpolat )
       mtparm->LoadDataToPair(CSD);
    else
       mtparm->LoadDataToLookup(CSD);
}

//=========================================================================================
// the same functions in ms_param.cpp



long int showMss = 1L;
/* test result GEM IPM calculation of equilibrium state in MULTI
long int TProfil::testMulti()
{
  MULTI *pmp = multi->GetPM();
  if( pmp->MK || pmp->PZ )
  {
   if( pa.p.PSM >= 2 )
   {
     fstream f_log("ipmlog.txt", ios::out|ios::app );
     f_log << "Warning " << pmp->stkey << ": " <<  pmp->errorCode << ":" << endl;
     f_log << pmp->errorBuf << endl;
   }
   if( showMss )
   {
           multi->addErrorMessage(" \nContinue?");
      switch( vfQuestion3(0, pmp->errorCode, pmp->errorBuf,
                           "&Yes", "&No", "&Yes to All" ))
       {
       case VF3_3:
           showMss=0l;
       case VF3_1:
           break;
       case VF3_2:
           Error(pmp->errorCode, pmp->errorBuf);
       }
   }

   return 1L;
  }

  return 0L	;
}
moved to TMulti*/

#ifdef USE_GEMS3K_SERVER
// Run process of calculate equilibria into the GEMSGUI shell
void  TProfil::CalculateEquilibriumGUI()
{
    TNodeGUI na(multi_internal);

#ifdef NO_ASYNC_SERVER
    zmq_req_client_t<TNodeGUI> zmqclient(na);
#else
    zmq_client_t<TNodeGUI> zmqclient(na);
#endif
    zmqclient.run_task();
}

#endif

namespace {
// FNV-1a hash of a block of memory.
void hashBytes( uint64_t& h, const void* data, size_t n )
{
    if( !data )
        return;
    const unsigned char* c = static_cast<const unsigned char*>( data );
    for( size_t i = 0; i < n; i++ ) { h ^= c[i]; h *= 1099511628211ULL; }
}

// True when the system has a phase with a kinetic model.
bool hasKineticPhases( const MULTI* pmp )
{
    if( !pmp->kMod )
        return false;
    for( long int k = 0; k < pmp->FI; k++ )
        if( pmp->kMod[k][0] != KM_UNDEF )
            return true;
    return false;
}

// True when TNode::unpackDataBr() would stop on the input: an element amount below pa_DB, or an upper
// species limit below the lower one. The direct solver call warns about these and continues.
bool inputRejectedByNode( const MULTI* pmp, double DB )
{
    for( long int i = 0; i < pmp->N - 1; i++ )
        if( pmp->B[i] < DB )
            return true;
    for( long int j = 0; j < pmp->L; j++ )
        if( pmp->DUL[j] < pmp->DLL[j] )
            return true;
    return false;
}

#ifdef USE_OPTIMA_SOLVER
// One kinetics step in AOP or SOP. If Optima fails, the step is solved again by the GEMS3K IPM solver
// from a cold start. The kinetic limits of the step are already set, so the retry does not run the kinetics again.
void optimaKineticsStep( TMulti* multi, long int& NumIterFIA, long int& NumIterIPM )
{
    MULTI* pmp = multi->GetPM();
    try {
        multi->CalculateEquilibriumStateOptima( NumIterFIA, NumIterIPM );
        return;
    }
    catch( TError& err ) {
        gui_logger->warn( "Kinetics time {}: Optima failed ({}); this step is solved by the GEMS3K IPM solver.",
                          pmp->kTau, err.mess );
    }
    // GEMSGUI_KINFAIL_EXPORT=<path>: writes the first failed step as GEMS3K files <path>-dat.lst (for a test case)
    static bool exported = false;
    if( const char* ex = std::getenv( "GEMSGUI_KINFAIL_EXPORT" ); ex && !exported )
    {
        exported = true;
        double Tai[4] = { pmp->TCc, pmp->TCc, 0., 0.1 }, Pai[4] = { pmp->Pc, pmp->Pc, 0., 0.1 };
        auto na = TNodeArrayGUI::create( 1, multi );
        na->setExportNode();
        na->MakeNodeStructuresOne( nullptr, true, Tai, Pai );
        ProcessProgressFunction msg = []( const std::string&, long ) { return false; };
        na->genGEMS3KInputFiles( std::string( ex ) + "-dat.lst", msg, 1, GEMS3KGenerator::default_type_f,
                                 false, false, false, false );
    }
    pmp->pNP = 0;
    multi->CalculateEquilibriumState( NumIterFIA, NumIterIPM, false );
}
#endif
} // namespace

// The node for GUI calculations: all species, phases and elements of the current system, one T and P point.
// It is made again when the names, the stoichiometry or the phase setup of the system change.
TNodeGUI& TProfil::calcNode()
{
    MULTI* pmp = multi_internal->GetPM();
    uint64_t h = 1469598103934665603ULL;
    const long int dims[4] = { pmp->N, pmp->L, pmp->FI, pmp->FIs };
    hashBytes( h, dims, sizeof( dims ) );
    hashBytes( h, pmp->SB, (size_t)pmp->N * ( MAXICNAME + MAXSYMB ) );
    hashBytes( h, pmp->SM, (size_t)pmp->L * MAXDCNAME );
    hashBytes( h, pmp->SF, (size_t)pmp->FI * ( MAXPHNAME + MAXSYMB ) );
    hashBytes( h, pmp->A, (size_t)pmp->N * pmp->L * sizeof( double ) );
    hashBytes( h, pmp->L1, (size_t)pmp->FI * sizeof( long int ) );
    hashBytes( h, pmp->PHC, (size_t)pmp->FI );
    hashBytes( h, pmp->DCC, (size_t)pmp->L );
    hashBytes( h, pmp->ICC, (size_t)pmp->N );
    hashBytes( h, pmp->MM, (size_t)pmp->L * sizeof( double ) );
    hashBytes( h, pmp->Awt, (size_t)pmp->N * sizeof( double ) );
    const char aalp = pmp->Aalp ? 1 : 0;
    hashBytes( h, &aalp, 1 );
    const std::string key = std::to_string( h );
    if( !calc_node || key != calc_node_key )
    {
        calc_node = std::make_shared<TNodeGUI>( multi_internal );
        // a wide T and P tolerance: the thermodynamic data come from the GUI, not from this one-point table
        double Tai[4] = { pmp->TCc, pmp->TCc, 0., 1e4 }, Pai[4] = { pmp->Pc, pmp->Pc, 0., 1e6 };
        const double fx = pmp->FX;
        if( pmp->FX == kTotalGibbsEnergyUnset )
            pmp->FX = 0.;   // packDataBr() would report the "not computed" marker as an error
        calc_node->MakeNodeStructures( nullptr, true, Tai, Pai );
        pmp->FX = fx;
        calc_node_key = key;
    }
    return *calc_node;
}

bool TProfil::processTitrates( bool runningOnly ) const
{
    const TProcess* proc = TProcess::pm;
    if( !proc || !proc->pep )
        return false;
    if( proc->pep->PsPro == S_OFF )   // only the inverse titration ('G') runs the golden section
        return false;
    if( runningOnly )
        return proc->pep->Istat == P_EXECUTE || proc->pep->Istat == P_MT_EXECUTE;
    return true;
}

std::string TProfil::targetsIgnoredReason( bool withProcess ) const
{
    const SYSTEM* syp = TSyst::sm->GetSY();
    if( syp->PpHtg != S_ON && syp->PEhtg != S_ON )
        return "";
    // a Process that titrates by itself does not use the targets in any solver mode: nothing to warn about
    if( withProcess && processTitrates( false ) )
        return "";
    if( !targetsUsedBy( solver_mode ) )
        return std::string( "The solver mode is " ) + solverModeName( solver_mode ) +
               ": the pH and Eh targets are not used. Select AOP, SOP, HOP or SHP in Calculate > Solver mode to use them.";
    return "";
}

// GEM IPM calculation of equilibrium state in MULTI
// without testing changes in the system
//
double TProfil::ComputeEquilibriumState( /*long int& NumPrecLoops,*/ long int& NumIterFIA, long int& NumIterIPM )
{
  TSysEq* STat = dynamic_cast<TSysEq *>(aMod[RT_SYSEQ].get());
  calcFinished = false;
  double calcTime = 0.;

  //multi->Access_GEM_IMP_init();
  //outMultiTxt( "Reaktoro_before.dump.txt"  );
#ifdef USE_GEMS3K_SERVER
  CalculateEquilibriumGUI( );
  calcTime = multi_internal->GetPM()->t_elap_sec;
#else
  // The Calculate menu sets the mode, cold or warm, for every calculation, also in a Process run.
  // A warm mode needs an earlier solution in memory.
  const SolverMode mode = solver_mode;
  MULTI* pmp = multi_internal->GetPM();
  // pH and Eh targets of the system recipe (Targets, tg_). Only the Optima modes use them: in AIA and
  // SIA, and in a build without Optima, a target is not used and a warning is written.
  const SYSTEM* syp = TSyst::sm->GetSY();
  const bool targetSet = syp->PpHtg == S_ON || syp->PEhtg == S_ON;
  // A Process that does its own titration (golden section) sets the amount itself: no target is used then.
  const bool processOwnTitration = processTitrates( true );
  const bool useTargets = targetSet && targetsUsedBy( mode ) && !processOwnTitration;
  const bool pHtarget = useTargets && syp->PpHtg == S_ON, Ehtarget = useTargets && syp->PEhtg == S_ON;
  {
      // warn once, and again only after the mode or the target setting changes
      static int lastWarned = -1;
      const std::string reason = processOwnTitration ? std::string() : targetsIgnoredReason( false );
      const int state = reason.empty() ? -1 : int( mode );
      if( state >= 0 && state != lastWarned )
          gui_logger->warn( "{}", reason );
      lastWarned = state;
  }
  struct TargetsOff {   // removes the targets after the calculation, also after an error
      TMulti* m;
      ~TargetsOff() { m->ClearControlConditions(); }
  } targetsOff{ multi_internal };
  multi_internal->ClearControlConditions();
  if( pHtarget || Ehtarget )
  {
#ifdef USE_OPTIMA_SOLVER
      // the target setters find the H, charge and H+ indices through the node
      TNodeGUI& tnode = calcNode();
      multi_internal->node1 = &tnode;
      try {
          if( pHtarget )
              multi_internal->SetControlCondition_pH( syp->pHtg );
          if( Ehtarget )
              multi_internal->SetControlCondition_Eh( syp->Ehtg );
      }
      catch( ... ) { multi_internal->node1 = nullptr; throw; }
      multi_internal->node1 = nullptr;
#endif
  }
  const bool warm = solverModeIsWarm( mode ) && pmp->pESU != 0;
#ifndef USE_OPTIMA_SOLVER
  if( mode != SolverMode::AIA && mode != SolverMode::SIA )
      Error( "Solver mode", "GEMS3K was built without the Optima solver (USE_OPTIMA_SOLVER): "
                            "choose AIA or SIA in Calculate > Solver mode." );
#endif
  // Kinetics, the step-by-step mode and an input that the node rejects use the direct solver call.
  const bool kinetics = pmp->kdT > 0. && hasKineticPhases( pmp );
  const bool rejected = inputRejectedByNode( pmp, multi_internal->base_param()->DB );
  if( kinetics || stepWise1 || rejected )
  {
      gui_logger->debug( "Direct solver call: kinetics {}, step-by-step {}, input rejected by the node {}",
                         kinetics, stepWise1, rejected );
      TMultiBase::BeginCalculation();
      pmp->pNP = warm ? 1 : 0;
      switch( mode )
      {
      case SolverMode::AIA: case SolverMode::SIA:
          multi_internal->CalculateEquilibriumState( NumIterFIA, NumIterIPM );
          break;
#ifdef USE_OPTIMA_SOLVER
      case SolverMode::AOP: case SolverMode::SOP:
          if( kinetics )
              optimaKineticsStep( multi_internal, NumIterFIA, NumIterIPM );
          else
              multi_internal->CalculateEquilibriumStateOptima( NumIterFIA, NumIterIPM );
          break;
      case SolverMode::HOP: case SolverMode::SHP:
          multi_internal->CalculateEquilibriumStateHOP( NumIterFIA, NumIterIPM, warm );
          break;
#endif
      default:
          break;
      }
      multi_internal->ReportAbsentPhaseIndexes();   // a real saturation index for a missing mixed phase
      calcTime = pmp->t_elap_sec;
  }
  else
  {
      // All other calculations run through TNode::GEM_run(), with the quick first try and the retries of GEMS3K.
      long int need = NEED_GEM_AIA;
      switch( mode )
      {
      case SolverMode::AIA: case SolverMode::SIA: need = warm ? NEED_GEM_SIA : NEED_GEM_AIA; break;
#ifdef USE_OPTIMA_SOLVER
      case SolverMode::AOP: case SolverMode::SOP: need = warm ? NEED_GEM_SOP : NEED_GEM_AOP; break;
      case SolverMode::HOP: case SolverMode::SHP: need = warm ? NEED_GEM_SHP : NEED_GEM_HOP; break;
#endif
      default: break;
      }
      TNodeGUI& node = calcNode();
      DATABR* dbr = node.pCNode();
      const DATACH* dch = node.pCSD();
      // GEM_run() writes these from the node; they are put back after, so the GUI keeps its own values
      const double TCc = pmp->TCc, Tc = pmp->Tc, Pc = pmp->Pc, VXc = pmp->VXc, kTau = pmp->kTau, kdT = pmp->kdT;
      // the inputs, from MULTI; the node keeps its own last answer, which a retry can start from
      dbr->TK = pmp->TCc + C_to_K;
      dbr->P = pmp->Pc * bar_to_Pa;
      dbr->Vs = pmp->VXc * 1e-6;
      for( long int i = 0; i < dch->nICb; i++ )
          dbr->bIC[i] = pmp->B[ dch->xic[i] ];
      for( long int j = 0; j < dch->nDCb; j++ )
      {
          dbr->dul[j] = pmp->DUL[ dch->xdc[j] ];
          dbr->dll[j] = pmp->DLL[ dch->xdc[j] ];
      }
      for( long int k = 0; k < dch->nPHb; k++ )
      {
          if( dch->nAalp > 0 )
              dbr->aPH[k] = pmp->Aalp[ dch->xph[k] ] * kg_to_g;
          dbr->omPH[k] = pmp->Falp[ dch->xph[k] ];
          dbr->xPH[k] = pmp->XF[ dch->xph[k] ];
      }
      node.setTargetsFromSystem();   // GEM_run() registers these again from the DBR
      if( !useTargets )
          dbr->pHtarget = dbr->Ehtarget = DOUBLE_EMPTY;   // not used in this mode
      dbr->Tm = 0.;   // no kinetic phases: the time step stays with the Process only
      dbr->dt = 0.;
      dbr->NodeStatusCH = need;
      // the retries must not ask questions; one question is asked at the end
      const long int showMss0 = showMss;
      showMss = 0L;
      const long int st = node.GEM_run( false );
      showMss = showMss0;
      node.SyncMultiWithDataBr();   // a retry can return an earlier attempt, kept only in the node
      pmp->TCc = TCc; pmp->Tc = Tc; pmp->Pc = Pc; pmp->VXc = VXc; pmp->kTau = kTau; pmp->kdT = kdT;
      long int numK2 = 0;
      calcTime = node.GEM_CalcTime( numK2, NumIterFIA, NumIterIPM );
      pmp->t_elap_sec = calcTime;
      gui_logger->debug( "Node solver call: mode {}, status {}, G {}", need, st, pmp->FX );
      // a cancelled call ends as a failure of any kind: it is a cancel, not an error to show
      const bool failed = st == T_ERROR_GEM || st == ERR_GEM_AIA || st == ERR_GEM_SIA || st == BAD_GEM_AIA || st == BAD_GEM_SIA
#ifdef USE_OPTIMA_SOLVER
          || st == ERR_GEM_AOP || st == ERR_GEM_SOP || st == ERR_GEM_HOP || st == ERR_GEM_SHP
          || st == BAD_GEM_AOP || st == BAD_GEM_SOP || st == BAD_GEM_HOP || st == BAD_GEM_SHP
#endif
          ;
      if( userCancel1 && failed )
          throw TMulti::UserCancelException();
      if( st == T_ERROR_GEM || st == ERR_GEM_AIA || st == ERR_GEM_SIA
#ifdef USE_OPTIMA_SOLVER
          || st == ERR_GEM_AOP || st == ERR_GEM_SOP || st == ERR_GEM_HOP || st == ERR_GEM_SHP
#endif
         )
          Error( "GEM IPM calculation", node.ipmLogError().empty() ? std::string( "no result" ) : node.ipmLogError() );
      if( st == BAD_GEM_AIA || st == BAD_GEM_SIA
#ifdef USE_OPTIMA_SOLVER
          || st == BAD_GEM_AOP || st == BAD_GEM_SOP || st == BAD_GEM_HOP || st == BAD_GEM_SHP
#endif
         )
          multi_internal->testMulti();   // logs the warning and asks "Continue?" as the direct call does
  }
#ifdef USE_OPTIMA_SOLVER
  // what each target needed: the H+ (pH) or electrons (Eh) added to or removed from the bulk composition
  for( const auto& c : multi_internal->GetControlConditions() )
  {
      const double added = -c.titrantAmount;   // H+ (pH) or electrons (Eh) added to the bulk composition
      if( c.targetMet )
          gui_logger->info( "Target {} = {}: reached {}; {} {} mol {}", c.name, c.target, c.achievedValue,
                            added >= 0. ? "added" : "removed", std::fabs( added ), c.name == "pH" ? "H+" : "electrons" );
      else
          gui_logger->warn( "Target {} = {} not reached: {} (tolerance {})", c.name, c.target, c.achievedValue, c.tolerance );
  }
#endif
#endif
  //outMultiTxt( "Reaktoro_after.dump.txt"  );

  calcFinished = true;
  STat->setCalcFlag( true );
  // STat->CellChanged(); // SD 28/11/2011 to protect MessageToSave()
  return calcTime;
}


// ------------------ End of m_param.cpp -----------------------