Skip to content

Commit e168360

Browse files
danv61danv61
authored
Merge pull request #1454 from virtualcell/dan-ss-rates-sanity
Langevin input file gets wrong particle location name
2 parents 2f3a714 + 0cefb78 commit e168360

5 files changed

Lines changed: 63 additions & 37 deletions

File tree

‎vcell-client/src/main/java/cbit/vcell/client/desktop/biomodel/ReactionRulePropertiesTableModel.java‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -241,6 +241,9 @@ public boolean isCellEditable(int rowIndex, int columnIndex) {
241241
* and the property that has changed.
242242
*/
243243
public void propertyChange(java.beans.PropertyChangeEvent evt) {
244+
// there is also an "entityChange" propertyChange event generated from ReactionRulePropertiesTreeModel (which is deprecated)
245+
// when a rule is selected in the upper panel (BioModelEditorReactionTableModel)
246+
// TODO: make sure we ignore it everywhere, should be phased out
244247
if (evt.getSource() == this && evt.getPropertyName().equals("reactionRule")) {
245248
ReactionRule oldValue = (ReactionRule)evt.getOldValue();
246249
if (oldValue != null){

‎vcell-client/src/main/java/cbit/vcell/mapping/gui/MolecularStructuresPropertiesPanel.java‎

Lines changed: 15 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -232,9 +232,23 @@ public void mouseClicked(MouseEvent e) {
232232
gbc.anchor = GridBagConstraints.WEST;
233233
optionsPanel.add(getZoomSmallerButton(), gbc);
234234

235+
// gbc = new GridBagConstraints();
236+
// gbc.gridx = 0;
237+
// gbc.gridy = 2;
238+
// gbc.anchor = GridBagConstraints.WEST;
239+
// gbc.insets = new Insets(4, 4, 4, 10);
240+
// optionsPanel.add(new JLabel("Reaction Radius"), gbc);
241+
//
242+
// gbc = new GridBagConstraints();
243+
// gbc.gridx = 1;
244+
// gbc.gridy = 2;
245+
// gbc.anchor = GridBagConstraints.WEST;
246+
// gbc.insets = new Insets(4, 4, 4, 10);
247+
// optionsPanel.add(new JLabel("2 nm"), gbc);
248+
235249
gbc = new GridBagConstraints();
236250
gbc.gridx = 0;
237-
gbc.gridy = 2;
251+
gbc.gridy = 2; // make this 3 if we show reaction radius here
238252
gbc.weightx = 1;
239253
gbc.weighty = 1; // fake cell used for filling all the vertical empty space
240254
gbc.anchor = GridBagConstraints.WEST;

‎vcell-client/src/main/java/org/vcell/model/springsalad/gui/MolecularStructuresPanel.java‎

Lines changed: 7 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -49,6 +49,7 @@
4949
import java.awt.event.FocusEvent;
5050
import java.awt.event.FocusListener;
5151
import java.beans.PropertyChangeListener;
52+
import java.text.DecimalFormat;
5253
import java.util.ArrayList;
5354
import java.util.Map;
5455

@@ -543,7 +544,12 @@ public Component getTableCellRendererComponent(JTable table, Object value, boole
543544
} else {
544545
setText(value + ""); // if it's too busy, just show the numbers
545546
}
546-
setToolTipText(value + " [nm]"); // we always show the units in the tooltip
547+
double minRadius = Math.max(0.5+dvalue, 1.5*dvalue);
548+
Double R = Math.min(minRadius, dvalue+2);
549+
double truncatedR = Math.floor(R * 100) / 100.0;
550+
DecimalFormat decimalFormat = new DecimalFormat("#.##");
551+
String formattedR = decimalFormat.format(truncatedR);
552+
setToolTipText(value + " [nm], (Reaction Radius: " + formattedR + " [nm])"); // we always show the units in the tooltip
547553
break;
548554
case COLUMN_DIFFUSION:
549555
if(cellWidth > 70) {

‎vcell-core/src/main/java/cbit/vcell/mapping/ReactionRuleSpec.java‎

Lines changed: 22 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -1201,12 +1201,14 @@ public void gatherIssues(IssueContext issueContext, List<Issue> issueList, React
12011201
if(scCandidate == null || scCandidate.getSpeciesPattern() == null) {
12021202
continue;
12031203
}
1204+
// we may have just one spCandidate wif the reaction is A + A -> A.A
12041205
SpeciesPattern spCandidate = scCandidate.getSpeciesPattern();
12051206
MolecularTypePattern mtpCandidate = spCandidate.getMolecularTypePatterns().get(0);
12061207
MolecularType mtCandidate = mtpCandidate.getMolecularType();
12071208
if(mtOursOne == mtCandidate) {
12081209
siteAttributesMapOne = scs.getSiteAttributesMap();
1209-
} else if(mtOursTwo == mtCandidate) {
1210+
}
1211+
if(mtOursTwo == mtCandidate) {
12101212
siteAttributesMapTwo = scs.getSiteAttributesMap();
12111213
}
12121214
}
@@ -1225,13 +1227,14 @@ public void gatherIssues(IssueContext issueContext, List<Issue> issueList, React
12251227
}
12261228
}
12271229
//
1228-
// TODO: check logic here, sometimes siteAttributesMapTwo is null (issue #977)
1229-
// temp fix, we check for non-null siteAttributesMapTwo
1230-
// hard to catch, it may be that the reactant or the product end up in an inconsistent state
1231-
// after using the reaction visual editor in the physiology (the state of a bonding site is possibly inconsistent)
1232-
//
1233-
if(sasOne != null && sasTwo != null && siteAttributesMapTwo != null) {
1234-
1230+
// with previous logic, siteAttributesMapTwo could have been null (issue #977)
1231+
// happened when the reaction was between 2 molecules of the same type A + A -> A.A
1232+
// the temp fixwas: we check for non-null siteAttributesMapTwo
1233+
// as from feb 24, 2025 this was fixed and siteAttributesMapTwo should never be null
1234+
if(siteAttributesMapTwo == null) {
1235+
throw new RuntimeException("Unexpected null value for siteAttributesMapTwo");
1236+
}
1237+
if(sasOne != null && sasTwo != null) {
12351238
for(Map.Entry<MolecularComponentPattern, SiteAttributesSpec> entry : siteAttributesMapTwo.entrySet()) {
12361239
MolecularComponentPattern mcpCandidate = entry.getKey();
12371240
if(MolecularComponentPattern.BondType.None != mcpCandidate.getBondType()) {
@@ -1250,12 +1253,12 @@ public void gatherIssues(IssueContext issueContext, List<Issue> issueList, React
12501253
issueList.add(new Issue(r, issueContext, IssueCategory.Identifiers, msg, tip, Issue.Severity.ERROR));
12511254
return;
12521255
}
1253-
// if(checkOnRate(sasOne, sasTwo) == false) { // rate doesn't check as acceptable
1254-
// String msg = "The simulation Kon is too large (I.e. exceeds the diffusion limited rate) for this reaction.";
1255-
// String tip = "Please consider reducing Kon or increasing the Radius or D of the participating Site Types.";
1256-
// issueList.add(new Issue(r, issueContext, IssueCategory.Identifiers, msg, tip, Issue.Severity.ERROR));
1257-
// return;
1258-
// }
1256+
if(checkOnRate(sasOne, sasTwo) == false) { // rate doesn't check as acceptable
1257+
String msg = "The forward rate Kf is too large (i.e. exceeds the diffusion limited rate) for this reaction rule.";
1258+
String tip = "Please consider reducing Kon or increasing the Radius or D of the participating Site Types.";
1259+
issueList.add(new Issue(r, issueContext, IssueCategory.Identifiers, msg, tip, Issue.Severity.WARNING));
1260+
return;
1261+
}
12591262
}
12601263
// binding reactions must be reversible
12611264
if(!reactionRule.isReversible()) {
@@ -1321,6 +1324,11 @@ public void gatherIssues(IssueContext issueContext, List<Issue> issueList, React
13211324
*/
13221325
public boolean checkOnRate(SiteAttributesSpec sasOne, SiteAttributesSpec sasTwo) {
13231326

1327+
// set of acceptable numbers (marginally) for A + A -> A.A are:
1328+
// Kon = 40 s-1uM-1
1329+
// site radius 1nm
1330+
// site diffusion rate 1 um2s-1
1331+
13241332
double R = sasOne.computeReactionRadius() + sasTwo.computeReactionRadius(); // nm
13251333
double D = sasOne.getDiffusionRate() + sasTwo.getDiffusionRate();
13261334

‎vcell-core/src/main/java/org/vcell/solver/langevin/LangevinLngvWriter.java‎

Lines changed: 16 additions & 21 deletions
Original file line numberDiff line numberDiff line change
@@ -5,6 +5,7 @@
55

66
import cbit.vcell.geometry.AnalyticSubVolume;
77
import cbit.vcell.geometry.SubVolume;
8+
import cbit.vcell.math.*;
89
import cbit.vcell.model.Structure;
910
import cbit.vcell.solver.*;
1011
import org.vcell.util.Pair;
@@ -18,28 +19,8 @@
1819
import cbit.vcell.mapping.SimulationContext;
1920
import cbit.vcell.mapping.SpeciesContextSpec;
2021
import cbit.vcell.mapping.SimulationContext.Application;
21-
import cbit.vcell.math.Action;
22-
import cbit.vcell.math.JumpProcessRateDefinition;
23-
import cbit.vcell.math.LangevinParticleJumpProcess;
24-
import cbit.vcell.math.LangevinParticleMolecularComponent;
25-
import cbit.vcell.math.LangevinParticleMolecularType;
26-
import cbit.vcell.math.MacroscopicRateConstant;
27-
import cbit.vcell.math.MathDescription;
28-
import cbit.vcell.math.MathException;
29-
import cbit.vcell.math.MathUtilities;
30-
import cbit.vcell.math.ParticleComponentStateDefinition;
31-
import cbit.vcell.math.ParticleComponentStatePattern;
32-
import cbit.vcell.math.ParticleJumpProcess;
33-
import cbit.vcell.math.ParticleMolecularComponent;
34-
import cbit.vcell.math.ParticleMolecularComponentPattern;
35-
import cbit.vcell.math.ParticleMolecularType;
36-
import cbit.vcell.math.ParticleMolecularTypePattern;
37-
import cbit.vcell.math.ParticleProperties;
3822
import cbit.vcell.math.ParticleProperties.ParticleInitialCondition;
3923
import cbit.vcell.math.ParticleProperties.ParticleInitialConditionCount;
40-
import cbit.vcell.math.ParticleSpeciesPattern;
41-
import cbit.vcell.math.SubDomain;
42-
import cbit.vcell.math.Variable;
4324
import cbit.vcell.mathmodel.MathModel;
4425
import cbit.vcell.messaging.server.SimulationTask;
4526
import cbit.vcell.parser.DivideByZeroException;
@@ -910,7 +891,21 @@ private static void writeSpeciesInfo(StringBuilder sb) {
910891
throw new RuntimeException("Initial concentration must be a number");
911892
}
912893

913-
sb.append("MOLECULE: \"" + lpmt.getName() + "\" " + subDomain.getName() +
894+
String locationName;
895+
if(subDomain instanceof MembraneSubDomain) {
896+
// we can't use the "Extracellular_Intracellular_membrane" math name of the subdomain,
897+
// the Langevin solver has no clue about what that might be
898+
// luckily we enforce that a membrane molecule must be anchored to exactly one membrane (named "Membrane")
899+
// alternatively we can check the Anchor site and get the membrane name "Membrane" from there
900+
if(lpmt.getAnchorList().size() != 1) {
901+
throw new RuntimeException("MembraneSubDomain particle must be anchored to one and only one Membrane");
902+
}
903+
// we make sure elsewhere that it's a membrane
904+
locationName = lpmt.getAnchorList().get(0);
905+
} else {
906+
locationName = subDomain.getName();
907+
}
908+
sb.append("MOLECULE: \"" + lpmt.getName() + "\" " + locationName +
914909
" Number " + scount +
915910
// number of site types and number of sites is the same for the vcell implementation of springsalad
916911
" Site_Types " + lpmt.getComponentList().size() + " Total" + "_Sites " + lpmt.getComponentList().size() +

0 commit comments

Comments
 (0)