genssi-developer / GenSSI

GenSSI 2.0 - Generating Series for testing Structural Identifiability
14 stars 3 forks source link

invalid index #20

Open craziecaesium opened 1 year ago

craziecaesium commented 1 year ago

Hi,

I am student working on a project regarding modelling of enzymatic networks. I am using GenSSI from Github: https://github.com/genssi-developer/GenSSI

I used the function genssiMain('rxnodes') and I got an error in the end, but I really don't understand how to resolve it.

It says:

Index in position 1 is invalid. Array indices must be positive integers or
logical values.

Error in genssiOrderTableau>displayRelevantParameters (line 291)
                    RJacparam_new(Mat_index(indX,j),:)=0;

Error in genssiOrderTableau (line 188)
                    displayRelevantParameters...

Error in genssiMain (line 137)
    [options,results] = genssiOrderTableau(model,results,RJacParam01,ECC,rParam,options);

 

I am using MATLAB R2022b if it is of any relevance.

Below is my model:

function model = rxnodes()

    % Specify stoichiometry matrix!
    R_matrix = [-1 0 0;1 -1 0;0 -1 0;0 -1 0;0 1 0;0 1 -1; 0 0 1];

    syms PRCA ACA NADPH CO2 NADP MMC SCA
    syms k_pco k_pca k_ccr k_aca k_nadph k_co2 k_mcm K_MMC
    syms pco ccr mcm

    % Parameters
    model.sym.p = [k_pco;k_pca;k_ccr;k_aca;k_nadph;k_co2;k_mcm;K_MMC];

    % State variables
    model.sym.x = [PRCA;ACA;NADPH;CO2;NADP;MMC;SCA];

    % Control vectors
    model.sym.g = [];

    % Specify the rate equation for the first reaction (pco):
    r1 = k_pco*pco*(PRCA/k_pca)/(1+PRCA/k_pca);

    % Specify the rate equation for the second reaction (ccr):
    r2_N = k_ccr*ccr*(ACA*NADPH*CO2)/(k_aca*k_nadph*k_co2);
    r2_D = (1+ACA/k_aca)*(1+NADPH/k_nadph)*(1+CO2/k_co2);
    r2 = r2_N/r2_D;

    % Specify the rate equation for the third reaction (mcm):
    r3 = k_mcm*mcm*(MMC/K_MMC)/(1+MMC/K_MMC);

    % Find the values of dC
    v = [r1;r2;r3];
    model.sym.xdot = R_matrix*v;

    % Initial conditions
    model.sym.x0 = [100;0;50;50;0;0;0];

    % Observables
    model.sym.y = [PRCA;NADPH;CO2;MMC;SCA];

end

Can this be resolved?

thomassligon commented 1 year ago

Hello Severo, I will look at this. Tom (Dr. Thomas S. Ligon) Tom@thomassligon.info

craziecaesium commented 1 year ago

Dear Dr. Ligon,

Please excuse me for the negligence, the error in line 291 mentioned above was able to be resolved using the files from branch issue19. But I would like to raise another error that I encountered: using the same model as above and applying default settings, the result showed the jacobian is rank deficient, thus I set Nder as 5, and the following error showed:

Index in position 1 exceeds array bounds. Index must not exceed 2.

Error in genssiOrderTableau>displayRelevantParameters (line 267)
        if sum_RJacParam01_nonzero_rows_t(indX,:)==length_Mat_index_t;

Error in genssiOrderTableau (line 188)
                    displayRelevantParameters...

Error in genssiMain (line 137)
    [options,results] = genssiOrderTableau(model,results,RJacParam01,ECC,rParam,options);

The same error occurs when I set Nder as 6, and I have no problem running through the examples that use higher Nder.

thomassligon commented 1 year ago

Hello Severo,

I will look at this as well, but I need to wait for a long-running calculation to finish on this PC.

craziecaesium commented 1 year ago

Thank you so much for the update, and the new feature is really neat!

thomassligon commented 1 year ago

Hello Severo, I don't fully understand this. You reported 2 errors, and the first was fixed by the fix created for Issue #19, which is now part of the master branch. However, you also reported a second error, which has not been fixed. Is the second error a problem for you?

craziecaesium commented 1 year ago

Sorry, I closed the issue because when I try to run the analysis at order 6 using the latest version, the error does not show up, but I just tried order 8, the same error (at line 267) came back again.