Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
23 changes: 23 additions & 0 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -81,6 +81,29 @@ Further examples can be found directly within this repository in the [Examples P
* How to run a [stochastic simulation](https://github.com/draeger-lab/SBSCL/blob/master/src/main/java/fern/Start.java)


## Using OptSolvX (LP) in SBSCL

### Run the OptSolvX demo

* In your IDE, run: `org.simulator.optsolvx.OptSolvXDemo`.
* The demo builds a tiny LP and solves it via OptSolvX (CommonsMath backend for now).
* Make sure the OptSolvX **jdk8** artifact is on the classpath (declared in SBSCL’s `pom.xml`).

### Enable debug logs

Create the adapter with `debug = true`:

```java
LPSolverAdapter solver =
new OptSolvXSolverAdapter(new CommonsMathSolver(), true);
```

### (Preview) SBML/FBC β†’ LP

An experimental entry point is available at `org.simulator.optsolvx.FbaToOptSolvX`.



### Comparison to Similar Libraries

To compare SBSCL to other simulation engines and to benchmark its predictions and results, a separate project, [SBSCL simulator comparison](https://github.com/matthiaskoenig/sbscl-simulator-comparison), is available.
Expand Down
8 changes: 8 additions & 0 deletions pom.xml
Original file line number Diff line number Diff line change
Expand Up @@ -471,6 +471,14 @@
</exclusion>
</exclusions>
</dependency>

<!-- OptSolvX -->
<dependency>
<groupId>org.optsolvx</groupId>
<artifactId>optsolvx</artifactId>
<version>0.1.0-SNAPSHOT</version>
<classifier>jdk8</classifier>
</dependency>
</dependencies>


Expand Down
194 changes: 194 additions & 0 deletions src/main/java/org/simulator/optsolvx/FbaToOptSolvX.java
Original file line number Diff line number Diff line change
@@ -0,0 +1,194 @@
package org.simulator.optsolvx;

import org.optsolvx.model.AbstractLPModel;
import org.optsolvx.model.Constraint;
import org.optsolvx.model.OptimizationDirection;
import org.sbml.jsbml.*;
import org.sbml.jsbml.ext.fbc.*;

import java.text.MessageFormat;
import java.util.LinkedHashMap;
import java.util.Map;
import java.util.logging.Logger;

/**
* SBML/FBC -> OptSolvX bridge (minimal functional version).
* Maps reactions to variables (with bounds), builds SΒ·v=0 constraints, and sets the active objective.
* NOTE: This intentionally ignores InitialAssignments and StoichiometryMath for now (TODO).
*/
public final class FbaToOptSolvX {
private static final Logger LOGGER = Logger.getLogger(FbaToOptSolvX.class.getName());

private FbaToOptSolvX() { /* no instances */ }

/** Build an OptSolvX LP model from SBML/FBC (v1 or v2). */
public static AbstractLPModel fromSBML(SBMLDocument doc) {
if (doc == null || !doc.isSetModel()) {
throw new IllegalArgumentException("SBMLDocument must contain a Model.");
}
final Model m = doc.getModel();

// Prepare LP model
final AbstractLPModel lp = new AbstractLPModel();

// Resolve FBC namespaces for level/version
final int level = doc.getLevel();
final int version = doc.getVersion();
final String fbcNSv1 = FBCConstants.getNamespaceURI(level, version, 1);
final String fbcNSv2 = FBCConstants.getNamespaceURI(level, version, 2);

// --- 1) Variables (reactions) with bounds ---------------------------------------------
// Default bounds if FBC has no values: 0 .. +inf (conventional for FBA)
final Map<String, Double> lb = new LinkedHashMap<>();
final Map<String, Double> ub = new LinkedHashMap<>();

for (Reaction r : m.getListOfReactions()) {
String rid = r.getId();
double lower = 0.0d;
double upper = Double.POSITIVE_INFINITY;

// FBC v2: bounds via FBCReactionPlugin lower/upper "Parameter" instances
if (r.isSetPlugin(fbcNSv2)) {
FBCReactionPlugin rp = (FBCReactionPlugin) r.getPlugin(fbcNSv2);
Parameter lpi = rp.getLowerFluxBoundInstance();
Parameter upi = rp.getUpperFluxBoundInstance();
if (lpi != null) lower = valueOf(lpi);
if (upi != null) upper = valueOf(upi);
}
// Store preliminary bounds; FBC v1 may override below via FluxBounds list
lb.put(rid, lower);
ub.put(rid, upper);
}

// FBC v1: FluxBounds at model level (override)
FBCModelPlugin mpV1 = (FBCModelPlugin) m.getPlugin(fbcNSv1);
if (mpV1 != null && mpV1.isSetListOfFluxBounds()) {
for (FluxBound fb : mpV1.getListOfFluxBounds()) {
String rid = fb.getReaction();
if (rid == null) continue;
switch (fb.getOperation()) { // Java 8 compatible switch
case GREATER_EQUAL:
lb.put(rid, fb.getValue());
break;
case LESS_EQUAL:
ub.put(rid, fb.getValue());
break;
case EQUAL:
lb.put(rid, fb.getValue());
ub.put(rid, fb.getValue());
break;
default:
LOGGER.warning("FBC v1: Unsupported FluxBound operation on " + rid);
}
}
}

// Add variables to LP
for (Reaction r : m.getListOfReactions()) {
String rid = r.getId();
double lower = nvl(lb.get(rid), 0.0d);
double upper = nvl(ub.get(rid), Double.POSITIVE_INFINITY);
lp.addVariable(rid, lower, upper);
}

// --- 2) Mass-balance constraints SΒ·v = 0 (ignore boundary species) ---------------------
for (Species s : m.getListOfSpecies()) {
boolean isBoundary = s.isSetBoundaryCondition() && s.getBoundaryCondition();
if (isBoundary) continue;

Map<String, Double> coeffs = new LinkedHashMap<>();
// For each reaction, accumulate stoichiometry of species s
for (Reaction r : m.getListOfReactions()) {
double sum = 0.0d;
// Reactants contribute negative stoichiometry
for (SpeciesReference sr : r.getListOfReactants()) {
if (s.getId().equals(sr.getSpecies())) {
sum += -stoich(sr);
}
}
// Products contribute positive stoichiometry
for (SpeciesReference sr : r.getListOfProducts()) {
if (s.getId().equals(sr.getSpecies())) {
sum += +stoich(sr);
}
}
if (sum != 0.0d) {
coeffs.put(r.getId(), sum);
}
}
if (!coeffs.isEmpty()) {
lp.addConstraint("mb_" + s.getId(), coeffs, Constraint.Relation.EQ, 0.0d);
} else {
// Optional: many species simply don't appear; keep quiet to avoid noisy logs
}
}

// --- 3) Objective (active FBC objective) -----------------------------------------------
Map<String, Double> obj = new LinkedHashMap<>();
OptimizationDirection dir = OptimizationDirection.MAXIMIZE; // sensible default

// Prefer FBC v2 at model level; fall back to v1
FBCModelPlugin mpV2 = (FBCModelPlugin) m.getPlugin(fbcNSv2);
FBCModelPlugin mp = (mpV2 != null) ? mpV2 : mpV1;

if (mp != null && mp.isSetActiveObjective()) {
Objective o = mp.getActiveObjectiveInstance();
if (o != null) {
// Direction
Objective.Type t = o.getType();
if (t == Objective.Type.MINIMIZE) {
dir = OptimizationDirection.MINIMIZE;
} else if (t == Objective.Type.MAXIMIZE) {
dir = OptimizationDirection.MAXIMIZE;
}
// Coefficients
for (FluxObjective fo : o.getListOfFluxObjectives()) {
String rid = fo.getReaction();
if (rid != null) obj.put(rid, fo.getCoefficient());
}
} else {
LOGGER.warning("No active FBC Objective instance set; using empty objective.");
}
} else {
LOGGER.warning("FBC Objective not set; using empty objective (objective=0).");
}

if (!obj.isEmpty()) {
lp.setObjective(obj, dir);
} else {
// Keep objective empty β†’ objective value will be 0.0; direction still recorded if needed
lp.setObjective(new LinkedHashMap<String, Double>(), dir);
}

// --- 4) Finalize & Log -----------------------------------------------------------------
lp.build();
LOGGER.info(MessageFormat.format(
"FbaToOptSolvX: built LP (vars={0}, cons={1}, objectiveVars={2}, dir={3})",
lp.getVariables().size(),
lp.getConstraints().size(),
obj.size(),
dir
));
return lp;
}

// --- Helpers -----------------------------------------------------------------------------

/** Get numeric value from a Parameter (defaults to 0.0 if unset). */
private static double valueOf(Parameter p) {
if (p == null) return 0.0d;
// Prefer explicit value; ignoring units & possible initial assignments for now (TODO)
return p.isSetValue() ? p.getValue() : 0.0d;
}

/** Null-safe value with default. */
private static double nvl(Double v, double def) {
return (v != null) ? v.doubleValue() : def;
}

/** Basic stoichiometry extraction (ignores StoichiometryMath & InitialAssignments for now). */
private static double stoich(SpeciesReference sr) {
if (sr == null) return 0.0d;
return sr.isSetStoichiometry() ? sr.getStoichiometry() : 1.0d;
}
}
50 changes: 50 additions & 0 deletions src/main/java/org/simulator/optsolvx/OptSolvXDemo.java
Original file line number Diff line number Diff line change
@@ -0,0 +1,50 @@
package org.simulator.optsolvx;

import org.optsolvx.model.AbstractLPModel;
import org.optsolvx.model.Constraint;
import org.optsolvx.model.OptimizationDirection;
import org.optsolvx.solver.CommonsMathSolver;
import org.optsolvx.solver.LPSolution;
import org.optsolvx.solver.LPSolverAdapter;

import java.util.HashMap;
import java.util.Map;

public class OptSolvXDemo {
public static void main(String[] args) {
// Build a tiny LP: maximize x + y
AbstractLPModel model = new AbstractLPModel();
model.addVariable("x", 0.0d, Double.POSITIVE_INFINITY);
model.addVariable("y", 0.0d, Double.POSITIVE_INFINITY);

// Objective: max x + y (Java 8 style map)
Map<String, Double> objective = new HashMap<>();
objective.put("x", 1.0d);
objective.put("y", 1.0d);
model.setObjective(objective, OptimizationDirection.MAXIMIZE);

// Constraints (Java 8 style maps)
Map<String, Double> c1 = new HashMap<>();
c1.put("x", 1.0d);
c1.put("y", 2.0d);
model.addConstraint("c1", c1, Constraint.Relation.LEQ, 4.0d);

Map<String, Double> c2 = new HashMap<>();
c2.put("x", 1.0d);
model.addConstraint("c2", c2, Constraint.Relation.LEQ, 3.0d);

// Finalize model
model.build();

// Solve via adapter, using CommonsMath as backend
LPSolverAdapter backend = new CommonsMathSolver();
LPSolverAdapter solver = new OptSolvXSolverAdapter(backend, /*debug=*/true);

LPSolution sol = solver.solve(model);

// Print result (expected optimum: x=3.0, y=0.5, objective=3.5)
System.out.println("Variables: " + sol.getVariableValues());
System.out.println("Objective: " + sol.getObjectiveValue());
System.out.println("Feasible: " + sol.isFeasible());
}
}
66 changes: 66 additions & 0 deletions src/main/java/org/simulator/optsolvx/OptSolvXSolverAdapter.java
Original file line number Diff line number Diff line change
@@ -0,0 +1,66 @@
package org.simulator.optsolvx;

import org.optsolvx.model.AbstractLPModel;
import org.optsolvx.solver.LPSolution;
import org.optsolvx.solver.LPSolverAdapter;

import java.util.Objects;
import java.text.MessageFormat;
import java.util.logging.Logger;

public final class OptSolvXSolverAdapter implements LPSolverAdapter {
private static final Logger LOGGER =
Logger.getLogger(OptSolvXSolverAdapter.class.getName());

private final LPSolverAdapter backend;
private final boolean debug;

/**
* @param backend concrete OptSolvX solver (e.g. CommonsMathSolver)
* @param debug enable verbose logging
*/
public OptSolvXSolverAdapter(LPSolverAdapter backend, boolean debug) {
this.backend = Objects.requireNonNull(backend, "backend must not be null");
this.debug = debug;
}

/** Convenience: debug disabled by default. */
public OptSolvXSolverAdapter(LPSolverAdapter backend) {
this(backend, false);
}

@Override
public LPSolution solve(AbstractLPModel model) {
if (model == null) {
throw new IllegalArgumentException("model must not be null");
}
if (!model.isBuilt()) {
throw new IllegalStateException("Model must be built() before solving.");
}

if (debug) {
LOGGER.info(MessageFormat.format(
"{0}: solving with {1} (vars={2}, cons={3})",
getClass().getSimpleName(),
backend.getClass().getSimpleName(),
model.getVariables().size(),
model.getConstraints().size()
));
}

long t0 = System.nanoTime();
LPSolution sol = backend.solve(model);
long dtMs = (System.nanoTime() - t0) / 1_000_000L;

if (debug) {
LOGGER.info(MessageFormat.format(
"{0}: result feasible={1}, objective={2}, time={3} ms",
getClass().getSimpleName(),
sol.isFeasible(),
sol.getObjectiveValue(),
dtMs
));
}
return sol;
}
}
Loading