Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Error :Jacobian matrix singular in NLEQ. #1233

Open
ZiyiiiLiu opened this issue Aug 27, 2024 · 11 comments
Open

Error :Jacobian matrix singular in NLEQ. #1233

ZiyiiiLiu opened this issue Aug 27, 2024 · 11 comments

Comments

@ZiyiiiLiu
Copy link

Hi,

I write a for loop to loop over many random parameter sets, for each parameter set it will produce 20 daughter parameter sets by small value mutation of each parameter, then I will plot the bifurcation diagram for each daughter parameter set using roadrunner.

I have used try_except loop for each parameter set to make the code junp the the next parameter set if it meets a error.
Also, I added
r = bf.model2te(model_sclc_2, ics=ics_1)
r.conservedMoietyAnalysis = True

Most parameter sets are running smoothly, but several of them will cause error:
�[35mError: Error :Jacobian matrix singular in NLEQ. Failed to converge to steady state. Check if Jacobian matrix is non-invertible or steady state solution does not exist.�[0m
/var/spool/slurm/spool/job2388885/slurm_script: line 14: 3672679 Aborted

The error will force the code to stop. Do you have any ways to make the code jump to the next parameter set but not stops at the error?

Thank you!
Zi

@luciansmith
Copy link

Those errors should be catchable in a try/except block. It's not exactly the same error, but this works for me:

import tellurium as te

r = te.loada("""
    -> S1; k1*S1
    k1 = 1
    S1 = 0
""")


for k1 in [-1, 0, 1]:
    r.resetAll()
    r.k1 = k1
    try:
        r.steadyState()
        print("Successful steady state found at", k1)
    except:
        print("Unsuccessful steady state at", k1)

gives the following output:

Successful steady state found at -1
Unsuccessful steady state at 0
Successful steady state found at 1

If this still isn't working for you, can you post the relevant bits of your code so we can see what's going on?

@hsauro
Copy link

hsauro commented Aug 27, 2024 via email

@ZiyiiiLiu
Copy link
Author

ZiyiiiLiu commented Aug 28, 2024

Yes, I used try-except but it does not work for this error: below is my part code: (pop mean parameter sets)
because the parameter sets and daughter sets are all produced randomly, so no example of that which will cause the error.

here the fitness_function includes the roadrunner part:

def fitness_function(pars): 
    m_sclc,r_sclc = simulate_sclc(pars) 
    autofunction(r_sclc)
   ......
    return fitness_score, bf_data, result, columns_after_l2_norm


for gen in range(num_generations):
        # Step 2: Evaluate the fitness of each individual in the population
        fitnesses = []
        well_pop=0
        individual_fitness_pairs = []  #for individual-fitness pairs
        for idx, pars in enumerate(pop):
            print('gen:',gen)
            print('run:', run)
            print(f'Index:{idx}')
            print(pars) #check why there was no fc>0 after the 6pm today.
            try:
             **fitness_score, bf_data, result,columns_after_l2_norm = fitness_function(pars)**
                fitnesses.append(fitness_score)
                print('ft score for sclc', fitness_score)
                if fitness_score>0:
                   plot_results(bf_data,gen, pars, fitness_score, run,columns_after_l2_norm)
                   if (gen == 0 and fitness_score > 1) or (gen > 0 and fitness_score > 1):
                       well_pop += 1 # Update well_pop count if ft>1
                   print('got well pop updated')
            except Exception as e:  
                fitness_score=0  
                fitnesses.append(0) #to make the length of fitnesses is the same with length of pop
                print('exception...................')
                continue    
            except:  
                print(f"An unexpected error occurred with parameters {pars}.")
                continue  # Skip to the next parameter set

@luciansmith
Copy link

So, 'simulate_sclc' is the function that calls roadrunner.steadyState()? And when you run the code, it doesn't print 'exception...' or 'An unexpected error...'?

Can you give us the full message you get from Python that ends with the "Jacobian matrix singular in NLEQ." error?

@ZiyiiiLiu
Copy link
Author

Hi Lucian,

The similate_sclc() includes the ode sets and m = r.simulate(0, 100, 1000) return m,r
(The full codes is the codes that I sent you for the memory issue, by the way, how about the memory issue? :-))

The autofunction(r_sclc) is below:
def autofunction(r):
auto = Plugin("tel_auto2000")
auto.setProperty("SBML", r.getCurrentSBML())
auto.setProperty("ScanDirection", 'Positive')
auto.setProperty("PrincipalContinuationParameter", "S1")
auto.setProperty("PreSimulation", "True")
auto.setProperty("PreSimulationDuration", 30)
auto.setProperty("RL0", 0)
auto.setProperty("RL1", 10)
auto.setProperty("NMX", 8000)
auto.setProperty("NPR", 3)
auto.setProperty("KeepTempFiles", True)
auto.setProperty("DS", 0.001)
auto.setProperty("DSMIN", 0.001)
auto.setProperty("DSMAX", 0.01)
auto.execute()
pts1 = auto.BifurcationPoints
if 1:
lbl1 = auto.BifurcationLabels
biData1 = auto.BifurcationData

I attached the full message here, hope this helps (there are some other errors but not affect the run to me) . Actually, I have some print() in the script, but they will only print results if the whole run complete, include the exception.........

myjob.o2388885_JacobianError.txt

@luciansmith
Copy link

I found an example model with the same failure to calculate the steady state, and can report that wrapping it in a try/except block works just fine:

import tellurium as te

r = te.loads(r"C:\Users\Lucian\Desktop\temp-biomodels\final\BIOMD0000000643\BIOMD0000000643_url.xml")

try:
    r.steadyState()
except Exception as e:
    print("Failure of steady state:", e)

Reports:

Failure of steady state: Jacobian matrix singular in NLEQ. Failed to converge to steady state. Check if Jacobian matrix is non-invertible or steady state solution does not exist.

I think you must be calling the steady state function somewhere that's not in that try/except block.

@ZiyiiiLiu
Copy link
Author

ZiyiiiLiu commented Aug 28, 2024

Thanks! I will try this later.
Could you please tell me the details about r.steadystate()?
I always use: (model_sclc_2 is the odes for the system)

    r = bf.model2te(model_sclc_2, ics=ics_1)
    r.conservedMoietyAnalysis = True 
    for k in pars:
        r[k] = pars[k]
    m = r.simulate(0, 100, 1000)

Are they the same thing?

@luciansmith
Copy link

r.simulate() and r.steadyState() are two very different functions. The first runs a simulation, and the second calculates the steady state (the point at which the system doesn't change). The error you're getting shouldn't arise from r.simulate() at all, which means as far as I can tell that you're calling r.steadyState() at a different point in your program that happens to be outside of your try/catch block?

The bit of script you post above looks fine to me: you set a bunch of parameter values, then simulate; that's pretty straightforward. But it shouldn't ever produce the error you're getting.

It's possible that you're calling a different function somewhere that itself calls 'steadyState' somehow? The error message complains about 'line 370' of 'genetic_algorithm_functions.py': what is that line?

(And sorry about the delay on the memory issue! It's still on our radar...)

@ZiyiiiLiu
Copy link
Author

Could you please share the r"C:\Users\Lucian\Desktop\temp-biomodels\final\BIOMD0000000643\BIOMD0000000643_url.xml to me if convenient? I want to know and try how to use the steadystate(). Thanks!

@luciansmith
Copy link

I mean, I picked it because it failed, not because it was a good example of running something to steady state ;-)

But! It's Biomodel 643:

https://www.ebi.ac.uk/biomodels/services/download/get-files/MODEL1707020000/4/BIOMD0000000643_url.xml

If you load it and call steadyState(), it fails with the same exception you're getting.

Here's another simpler example, courtesy of Herbert:

import tellurium as te
r = te.loada("""
     A -> B; k1
     B -> ; k2
     A = 10
     k1 = 0.1; k2 = 0.2
""")

try:
    r.steadyState()
except Exception as e:
    print("Failure of steady state:", e)

@ZiyiiiLiu
Copy link
Author

ZiyiiiLiu commented Sep 2, 2024

I can print the error using your code, but my own code can't print the error and can't jump the error when going through auto.execute(), so the process will stop at here. (This is running in linux server)

But I tried the same code with same error in windows system, the whole code running smoothly, which can print the error and jump to another parameter set.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants