/*
 * postbuckling_shell_arclength.java
 */

import com.comsol.model.*;
import com.comsol.model.util.*;

/** Model exported on May 15 2026, 08:13 by COMSOL 6.4.0.421. */
public class postbuckling_shell_arclength {

  public static Model run() {
    Model model = ModelUtil.create("Model");

//    In this example, you will start from an existing model in the Structural Mechanics Module.
//    From the File menu, choose Application Libraries.
//    In the Application Libraries window, select Structural Mechanics Module > Verification Examples > postbuckling_shell in the tree.
//    Click Open.

    model.component().create("comp1", true);

    model.component("comp1").geom().create("geom1", 3);
    model.component("comp1").geom("geom1").geomRep("comsol");

    model.component("comp1").mesh().create("mesh1");
    model.component("comp1").mesh("mesh1").contribute("geom/detail", true);

    model.component("comp1").physics().create("shell", "Shell", "geom1");

    model.study().create("std1");
    model.study("std1").create("stat", "Stationary");

    model.param().set("R", "2540[mm]");
    model.param().descr("R", "Panel radius");
    model.param().set("L", "254[mm]");
    model.param().descr("L", "Panel length");
    model.param().set("thic", "6.35[mm]");
    model.param().descr("thic", "Panel thickness");
    model.param().set("theta", "0.1[rad]");
    model.param().descr("theta", "Panel section angle");
    model.param().set("E0", "3.103[GPa]");
    model.param().descr("E0", "Young's modulus");
    model.param().set("nu0", "0.3");
    model.param().descr("nu0", "Poisson's ratio");
    model.param().set("disp", "0");
    model.param().descr("disp", "Displacement parameter");

    model.component("comp1").geom("geom1").create("wp1", "WorkPlane");
    model.component("comp1").geom("geom1").feature("wp1").set("unite", true);
    model.component("comp1").geom("geom1").feature("wp1").set("quickplane", "xz");
    model.component("comp1").geom("geom1").feature("wp1").geom().create("ls1", "LineSegment");
    model.component("comp1").geom("geom1").feature("wp1").geom().feature("ls1").set("specify1", "coord");
    model.component("comp1").geom("geom1").feature("wp1").geom().feature("ls1").set("specify2", "coord");
    model.component("comp1").geom("geom1").feature("wp1").geom().feature("ls1")
         .set("coord1", new String[]{"0", "R"});
    model.component("comp1").geom("geom1").feature("wp1").geom().feature("ls1")
         .set("coord2", new String[]{"L", "R"});
    model.component("comp1").geom("geom1").feature("wp1").geom().run("ls1");
    model.component("comp1").geom("geom1").run("wp1");
    model.component("comp1").geom("geom1").feature().create("rev1", "Revolve");
    model.component("comp1").geom("geom1").feature("rev1").set("angtype", "specang");
    model.component("comp1").geom("geom1").feature("rev1").set("angle2", "theta");
    model.component("comp1").geom("geom1").feature("rev1").set("axis", new double[]{1, 0});
    model.component("comp1").geom("geom1").run("rev1");

    model.component("comp1").cpl().create("aveop1", "Average");

    model.component("comp1").geom("geom1").run();

    model.component("comp1").cpl("aveop1").set("axisym", true);
    model.component("comp1").cpl("aveop1").selection().geom("geom1", 2);
    model.component("comp1").cpl("aveop1").selection().set(1);
    model.component("comp1").cpl().create("intop1", "Integration");
    model.component("comp1").cpl("intop1").set("axisym", true);
    model.component("comp1").cpl("intop1").selection().geom("geom1", 0);
    model.component("comp1").cpl("intop1").selection().set(4);

    model.component("comp1").variable().create("var1");
    model.component("comp1").variable("var1").set("w_center", "-intop1(w)");
    model.component("comp1").variable("var1").descr("w_center", "Vertical displacement at shell center");

    model.component("comp1").physics("shell").feature("to1").set("d", "thic");
    model.component("comp1").physics("shell").create("sym1", "SymmetrySolid1", 1);
    model.component("comp1").physics("shell").feature("sym1").selection().set(3, 4);
    model.component("comp1").physics("shell").create("pin1", "Pinned", 1);
    model.component("comp1").physics("shell").feature("pin1").selection().set(2);
    model.component("comp1").physics("shell").create("pl1", "PointLoad", 0);
    model.component("comp1").physics("shell").feature("pl1").selection().set(4);
    model.component("comp1").physics("shell").feature("pl1").set("forcePoint", new String[]{"0", "0", "-P/4"});
    model.component("comp1").physics("shell").create("ge1", "GlobalEquations", -1);
    model.component("comp1").physics("shell").feature("ge1").setIndex("name", "P", 0, 0);
    model.component("comp1").physics("shell").feature("ge1").setIndex("equation", "aveop1(-w)-disp", 0, 0);
    model.component("comp1").physics("shell").feature("ge1").setIndex("description", "Force at shell center", 0, 0);
    model.component("comp1").physics("shell").feature("ge1").set("DependentVariableQuantity", "force");
    model.component("comp1").physics("shell").feature("ge1").set("SourceTermQuantity", "displacement");

    model.component("comp1").material().create("mat1", "Common");
    model.component("comp1").material("mat1").propertyGroup()
         .create("Enu", "Enu", "Young's_modulus_and_Poisson's_ratio");
    model.component("comp1").material("mat1").propertyGroup("Enu").set("E", new String[]{"E0"});
    model.component("comp1").material("mat1").propertyGroup("Enu").set("nu", new String[]{"nu0"});
    model.component("comp1").material("mat1").propertyGroup("def").set("density", new String[]{"0"});

    model.component("comp1").mesh("mesh1").create("map1", "Map");
    model.component("comp1").mesh("mesh1").feature("map1").selection().set(1);
    model.component("comp1").mesh("mesh1").feature("map1").create("dis1", "Distribution");
    model.component("comp1").mesh("mesh1").feature("map1").feature("dis1").selection().set(1, 2);
    model.component("comp1").mesh("mesh1").feature("map1").feature("dis1").set("numelem", 10);
    model.component("comp1").mesh("mesh1").run("map1");

    model.study("std1").label("Postbuckling Study");
    model.study("std1").feature("stat").set("useparam", true);
    model.study("std1").feature("stat").setIndex("pname", "R", 0);
    model.study("std1").feature("stat").setIndex("plistarr", "", 0);
    model.study("std1").feature("stat").setIndex("punit", "m", 0);
    model.study("std1").feature("stat").setIndex("pname", "R", 0);
    model.study("std1").feature("stat").setIndex("plistarr", "", 0);
    model.study("std1").feature("stat").setIndex("punit", "m", 0);
    model.study("std1").feature("stat").setIndex("pname", "disp", 0);
    model.study("std1").feature("stat").setIndex("plistarr", "range(0,2e-4,1)", 0);
    model.study("std1").feature("stat").set("geometricNonlinearity", true);
    model.study("std1").showAutoSequences("all");

    model.sol("sol1").feature("s1").feature("p1").create("st1", "StopCondition");
    model.sol("sol1").feature("s1").feature("p1").feature("st1").setIndex("stopcondarr", "", 0);
    model.sol("sol1").feature("s1").feature("p1").feature("st1").setIndex("stopcondterminateon", "true", 0);
    model.sol("sol1").feature("s1").feature("p1").feature("st1").setIndex("stopcondActive", true, 0);
    model.sol("sol1").feature("s1").feature("p1").feature("st1").setIndex("stopconddesc", "Stop expression 1", 0);
    model.sol("sol1").feature("s1").feature("p1").feature("st1").setIndex("stopcondarr", "", 0);
    model.sol("sol1").feature("s1").feature("p1").feature("st1").setIndex("stopcondterminateon", "true", 0);
    model.sol("sol1").feature("s1").feature("p1").feature("st1").setIndex("stopcondActive", true, 0);
    model.sol("sol1").feature("s1").feature("p1").feature("st1").setIndex("stopconddesc", "Stop expression 1", 0);
    model.sol("sol1").feature("s1").feature("p1").feature("st1").setIndex("stopcondarr", "comp1.w_center>0.035", 0);
    model.sol("sol1").feature("s1").feature("p1").feature("st1").set("storestopcondsol", "stepbefore");
    model.sol("sol1").feature("s1").feature("p1").feature("st1").set("stopcondwarn", false);
    model.sol("sol1").feature("s1").set("reacf", false);
    model.sol("sol1").runAll();

    model.result().dataset().create("dset1shellshl", "Shell");
    model.result().dataset("dset1shellshl").set("data", "dset1");
    model.result().dataset("dset1shellshl").setIndex("topconst", "1", 3, 1);
    model.result().dataset("dset1shellshl").setIndex("bottomconst", "-1", 3, 1);
    model.result().dataset("dset1shellshl").setIndex("orientationexpr", "shell.nlX", 0);
    model.result().dataset("dset1shellshl").setIndex("displacementexpr", "arx", 0);
    model.result().dataset("dset1shellshl").setIndex("orientationexpr", "shell.nlY", 1);
    model.result().dataset("dset1shellshl").setIndex("displacementexpr", "ary", 1);
    model.result().dataset("dset1shellshl").setIndex("orientationexpr", "shell.nlZ", 2);
    model.result().dataset("dset1shellshl").setIndex("displacementexpr", "arz", 2);
    model.result().dataset("dset1shellshl").set("distanceexpr", "shell.z_pos");
    model.result().dataset("dset1shellshl").set("seplevels", false);
    model.result().dataset("dset1shellshl").set("resolution", 2);
    model.result().dataset("dset1shellshl").set("areascalefactor", "shell.ASF");
    model.result().dataset("dset1shellshl").set("linescalefactor", "shell.LSF");
    model.result().create("pg1", "PlotGroup3D");
    model.result("pg1").set("data", "dset1shellshl");
    model.result("pg1").setIndex("looplevel", 92, 0);
    model.result("pg1").label("Stress (shell)");
    model.result("pg1").set("showlegends", true);
    model.result("pg1").set("frametype", "spatial");
    model.result("pg1").create("surf1", "Surface");
    model.result("pg1").feature("surf1").set("expr", new String[]{"shell.misesGp"});
    model.result("pg1").feature("surf1").set("threshold", "manual");
    model.result("pg1").feature("surf1").set("thresholdvalue", 0.2);
    model.result("pg1").feature("surf1").set("colortable", "Rainbow");
    model.result("pg1").feature("surf1").set("colortabletrans", "none");
    model.result("pg1").feature("surf1").set("colorscalemode", "linear");
    model.result("pg1").feature("surf1").set("descr", "von Mises stress");
    model.result("pg1").feature("surf1").set("colortable", "Prism");
    model.result("pg1").feature("surf1").create("def", "Deform");
    model.result("pg1").feature("surf1").feature("def").set("expr", new String[]{"shell.u", "shell.v", "shell.w"});
    model.result("pg1").feature("surf1").feature("def").set("scaleactive", true);
    model.result("pg1").feature("surf1").feature("def").set("scale", "1");
    model.result("pg1").run();
    model.result().evaluationGroup().create("eg1", "EvaluationGroup");
    model.result().evaluationGroup("eg1").label("Evaluation Group: Force vs. Displacement");
    model.result().evaluationGroup("eg1").create("pev1", "EvalPoint");
    model.result().evaluationGroup("eg1").feature("pev1").selection().set(4);
    model.result().evaluationGroup("eg1").feature("pev1").setIndex("expr", "w_center", 0);
    model.result().evaluationGroup("eg1").feature("pev1").setIndex("expr", "P", 1);
    model.result().evaluationGroup("eg1").set("includeparameters", false);
    model.result().evaluationGroup("eg1").run();
    model.result().create("pg2", "PlotGroup1D");
    model.result("pg2").run();
    model.result("pg2").label("Force vs. Displacement");
    model.result("pg2").set("xlabelactive", true);
    model.result("pg2").set("xlabel", "Vertical displacement at shell center (m)");
    model.result("pg2").set("ylabelactive", true);
    model.result("pg2").set("ylabel", "Force at shell center (N)");
    model.result("pg2").create("tblp1", "Table");
    model.result("pg2").feature("tblp1").set("markerpos", "datapoints");
    model.result("pg2").feature("tblp1").set("linewidth", "preference");
    model.result("pg2").feature("tblp1").set("source", "evaluationgroup");
    model.result("pg2").run();
    model.result("pg2").run();

    model.title("Postbuckling Analysis of a Hinged Cylindrical Shell");

    model
         .description("This example shows how to trace a postbuckling path where neither load nor displacement is increasing monotonously. The results are compared to published values.");

    model.label("postbuckling_shell.mph");

    model.result("pg2").run();

//    Remove the <l>Postbuckling Study</l>. Start by copying the evaluation group results to a table for future comparison.
//    In the Model Builder window, expand the Component 1 (comp1) node, then click Results > Evaluation Group: Force vs. Displacement.
//    In the Evaluation Group: Force vs. Displacement toolbar, click Copy to Table.

    model.result().table().create("tbl1", "Table");
    model.result().evaluationGroup("eg1").copyToTable("tbl1");

//    In the Settings window for Table, type Table: Reference Result in the Label text field.

    model.result().table("tbl1").label("Table: Reference Result");

//    In the Model Builder window, right-click Postbuckling Study and choose Delete.

    model.study().remove("std1");

//    Modify the parameter list.
//    In the Model Builder window, under Global Definitions, click Parameters 1.
//    In the Settings window for Parameters, locate the Parameters section.
//    In the table, enter the following settings:

    model.param().rename("disp", "p");
    model.param().descr("p", "Dimensionless load");
    model.param().set("p_scaleFactor", "1E5[N]");
    model.param().descr("p_scaleFactor", "Load scale factor");
    model.param().set("P", "p*p_scaleFactor");
    model.param().descr("P", "Total load");

//    In the Model Builder window, expand the Component 1 (comp1) > Definitions node.
//    Right-click Component 1 (comp1) > Definitions > Average 1 (aveop1) and choose Delete.

    model.component("comp1").cpl().remove("aveop1");

//    Delete the <l>Global Equation</l> node from the physics.
//    In the Model Builder window, expand the Component 1 (comp1) > Shell (shell) node.
//    Right-click Component 1 (comp1) > Shell (shell) > Global Equations 1 (ODE1) and choose Delete.

    model.component("comp1").physics("shell").feature().remove("ge1");

//    Create model methods to run the postbuckling study with the arc length method. Note that the method editor is only available in the Windows\[\textregistered\] version of the COMSOL Desktop.
//    In the Home toolbar, click Application Builder.
//    In the Home toolbar, click More Libraries and choose Utility Class.
//    Right-click util1 and choose Rename.
//    In the Rename Utility Class dialog, type codeUtil in the New name text field.
//    Click OK.
//    Right-click codeUtil and choose Edit.
//    Copy the following code into the <l>codeUtil</l> window:
//    //This utility creates the requested methods&#x2424;&#x2424;//Method to remove autocreated model nodes&#x2424;public static void removeAutoCreatedNodes() {&#x2424; if (contains(model.study().tags(), &quot;std_in&quot;))&#x2424; model.study().remove(&quot;std_in&quot;);&#x2424; if (contains(model.study().tags(), &quot;std_c&quot;))&#x2424; model.study().remove(&quot;std_c&quot;);&#x2424; if (contains(model.study().tags(), &quot;std_p&quot;))&#x2424; model.study().remove(&quot;std_p&quot;);&#x2424; if (contains(model.sol().tags(), &quot;sol_1&quot;))&#x2424; model.sol().remove(&quot;sol_1&quot;);&#x2424; if (contains(model.sol().tags(), &quot;sol_2&quot;))&#x2424; model.sol().remove(&quot;sol_2&quot;);&#x2424; if (contains(model.sol().tags(), &quot;sol_3&quot;))&#x2424; model.sol().remove(&quot;sol_3&quot;);&#x2424; if (contains(model.sol().tags(), &quot;sol_4&quot;))&#x2424; model.sol().remove(&quot;sol_4&quot;);&#x2424; if (contains(model.result().dataset().tags(), &quot;dset_c&quot;))&#x2424; model.result().dataset().remove(&quot;dset_c&quot;);&#x2424; if (contains(model.result().evaluationGroup().tags(), &quot;eg_c&quot;))&#x2424; model.result().evaluationGroup().remove(&quot;eg_c&quot;);&#x2424;}&#x2424;&#x2424;//Method to create a placeholder study for the initial guess&#x2424;public static void createPlaceholderStudyForInitialGuess() {&#x2424; Study study = model.study().create(&quot;std_in&quot;);&#x2424; String stdLabel = &quot;Placeholder Study for Initial Guess&quot;;&#x2424; boolean isLabelExists = false;&#x2424; for (Study std : model.study()) {&#x2424; if (std.label().equals(stdLabel))&#x2424; isLabelExists = true;&#x2424; }&#x2424; if (!isLabelExists)&#x2424; study.label(stdLabel);&#x2424; study.create(&quot;stat_in&quot;, &quot;Stationary&quot;);&#x2424; study.feature(&quot;stat_in&quot;).set(&quot;geometricNonlinearity&quot;, true);&#x2424; SolverSequence solin = model.sol().create(&quot;sol_in&quot;);&#x2424; solin.createAutoSequence(&quot;std_in&quot;);&#x2424; solin.attach(&quot;std_in&quot;);&#x2424;}&#x2424;&#x2424;//Method to create a study for the postbuckling analysis&#x2424;public static void createPostbucklingStudy() {&#x2424; Study study = model.study().create(&quot;std_c&quot;);&#x2424; String stdLabel = &quot;Postbuckling Study&quot;;&#x2424; boolean isLabelExists = false;&#x2424; for (Study std : model.study()) {&#x2424; if (std.label().equals(stdLabel))&#x2424; isLabelExists = true;&#x2424; }&#x2424; if (!isLabelExists)&#x2424; study.label(stdLabel);&#x2424; study.create(&quot;stat_c&quot;, &quot;Stationary&quot;);&#x2424; study.feature(&quot;stat_c&quot;).set(&quot;geometricNonlinearity&quot;, true);&#x2424; SolverSequence solc = model.sol().create(&quot;sol_c&quot;);&#x2424; solc.createAutoSequence(&quot;std_c&quot;);&#x2424; solc.attach(&quot;std_c&quot;);&#x2424;}&#x2424;&#x2424;//Method to create a placeholder study for the parametric solution&#x2424;public static void createPlaceholderStudyForParametricSolution(String p) {&#x2424; Study study = model.study().create(&quot;std_p&quot;);&#x2424; String stdLabel = &quot;Placeholder Study for Parametric Solution&quot;;&#x2424; boolean isLabelExists = false;&#x2424; for (Study std : model.study()) {&#x2424; if (std.label().equals(stdLabel))&#x2424; isLabelExists = true;&#x2424; }&#x2424; if (!isLabelExists)&#x2424; study.label(stdLabel);&#x2424; study.create(&quot;stat_p&quot;, &quot;Stationary&quot;);&#x2424; SolverSequence solp = model.sol().create(&quot;sol_p&quot;);&#x2424; solp.createAutoSequence(&quot;std_p&quot;);&#x2424; solp.attach(&quot;std_p&quot;);&#x2424; StudyFeature studyFeature = model.study(&quot;std_p&quot;).feature(&quot;stat_p&quot;);&#x2424; studyFeature.set(&quot;geometricNonlinearity&quot;, true);&#x2424; studyFeature.set(&quot;useparam&quot;, true);&#x2424; studyFeature.setIndex(&quot;pname&quot;, p, 0);&#x2424; studyFeature.setIndex(&quot;plistarr&quot;, &quot;1 2&quot;, 0);&#x2424; studyFeature.setIndex(&quot;punit&quot;, &quot;N&quot;, 0);&#x2424;}&#x2424;&#x2424;//Method to create an evaluation group&#x2424;public static void createEvaluationGroup(String globalVariable) {&#x2424; DatasetFeature dsetc = model.result().dataset().create(&quot;dset_c&quot;, &quot;Solution&quot;);&#x2424; dsetc.set(&quot;solution&quot;, &quot;sol_c&quot;);&#x2424; EvaluationGroupFeature eg = model.result().evaluationGroup().create(&quot;eg_c&quot;, &quot;EvaluationGroup&quot;);&#x2424; eg.set(&quot;data&quot;, &quot;dset_c&quot;);&#x2424; eg.create(&quot;gev1&quot;, &quot;EvalGlobal&quot;);&#x2424; eg.feature(&quot;gev1&quot;).setIndex(&quot;expr&quot;, globalVariable, 0);&#x2424; eg.set(&quot;includeparameters&quot;, false);&#x2424;}&#x2424;&#x2424;//Method for creating a stop condition in the solver&#x2424;public static boolean stopCondition(double globalVariableMin, double globalVariableMax) {&#x2424; boolean stop = false;&#x2424; EvaluationGroupFeature eg = model.result().evaluationGroup(&quot;eg_c&quot;);&#x2424; eg.run();&#x2424; double[][] globalVariable = eg.getReal();&#x2424; if (globalVariable[0][0] &gt; globalVariableMax || globalVariable[0][0] &lt; globalVariableMin)&#x2424; stop = true;&#x2424; return stop;&#x2424;}&#x2424;
//    In the Application Builder window, right-click Methods and choose New Method.
//    In the New Method dialog, type postBucklingWithArclength in the Name text field.
//    Click OK.
//    In the Application Builder window, under Methods, click postBucklingWithArclength.
//    Copy the following code into the <l>postBucklingWithArclength</l> window:
//    clearDebugLog();&#x2424;long startTime = timeStamp();&#x2424;if (!contains(model.param().varnames(), p))&#x2424; error(&quot;The load parameter with the given name does not exist in the parameter list of the model, check the model.&quot;);&#x2424;&#x2424;if (p.equals(&quot;&quot;) || deltap == 0 || numberOfIter == 0 || (isStopConditionGiven &amp;&amp; globalVariable.equals(&quot;&quot;)))&#x2424; error(&quot;The inputs assigned to the model method are either empty or zero, check the model.&quot;);&#x2424;&#x2424;if (isStopConditionGiven &amp;&amp; globalVariableMin == globalVariableMax)&#x2424; error(&quot;The lower and upper limits of the global variable in stop condition should not be the same, check the model.&quot;);&#x2424;&#x2424;model.param().set(p, &quot;0&quot;);&#x2424;//Remove the autocreated model nodes&#x2424;codeUtil.removeAutoCreatedNodes();&#x2424;//Create a placeholder study to store the intial guess&#x2424;codeUtil.createPlaceholderStudyForInitialGuess();&#x2424;//Create a postbuckling study&#x2424;codeUtil.createPostbucklingStudy();&#x2424;//Create a placeholder study to store the parametric solutions&#x2424;codeUtil.createPlaceholderStudyForParametricSolution(p);&#x2424;//Create an evaluation group for the evaluation of the stop condition&#x2424;if (isStopConditionGiven)&#x2424; codeUtil.createEvaluationGroup(globalVariable);&#x2424;&#x2424;//First iteration&#x2424;model.sol(&quot;sol_c&quot;).runAll();&#x2424;model.sol(&quot;sol_p&quot;).setU(1, model.sol(&quot;sol_c&quot;).getU());&#x2424;model.sol(&quot;sol_c&quot;).copySolution(&quot;sol_1&quot;);&#x2424;//Second iteration&#x2424;model.param().set(p, toString(deltap));&#x2424;model.sol(&quot;sol_c&quot;).runAll();&#x2424;model.sol(&quot;sol_p&quot;).setU(2, model.sol(&quot;sol_c&quot;).getU());&#x2424;model.sol(&quot;sol_c&quot;).copySolution(&quot;sol_2&quot;);&#x2424;//Third iteration&#x2424;model.param().set(p, toString(2*deltap));&#x2424;model.sol(&quot;sol_c&quot;).runAll();&#x2424;model.sol(&quot;sol_p&quot;).setU(3, model.sol(&quot;sol_c&quot;).getU());&#x2424;model.sol(&quot;sol_c&quot;).copySolution(&quot;sol_3&quot;);&#x2424;//Fourth iteration&#x2424;model.param().set(p, toString(3*deltap));&#x2424;model.sol(&quot;sol_c&quot;).runAll();&#x2424;model.sol(&quot;sol_p&quot;).setU(4, model.sol(&quot;sol_c&quot;).getU());&#x2424;model.sol(&quot;sol_c&quot;).copySolution(&quot;sol_4&quot;);&#x2424;&#x2424;//Next iterations by the arclength method&#x2424;double[] parametricFull = new double[4+numberOfIter];&#x2424;parametricFull[0] = 0;&#x2424;parametricFull[1] = 1*deltap;&#x2424;parametricFull[2] = 2*deltap;&#x2424;parametricFull[3] = 3*deltap;&#x2424;int pfinal = 4+numberOfIter;&#x2424;model.study(&quot;std_c&quot;).feature(&quot;stat_c&quot;).set(&quot;useinitsol&quot;, true);&#x2424;model.study(&quot;std_c&quot;).feature(&quot;stat_c&quot;).set(&quot;initmethod&quot;, &quot;sol&quot;);&#x2424;model.study(&quot;std_c&quot;).feature(&quot;stat_c&quot;).set(&quot;initstudy&quot;, &quot;std_in&quot;);&#x2424;model.study(&quot;std_c&quot;).feature(&quot;stat_c&quot;).set(&quot;initsol&quot;, &quot;sol_in&quot;);&#x2424;int sizeOfSolution = pfinal;&#x2424;&#x2424;//Arclength method&#x2424;for (int ii = 4; ii &lt; pfinal; ii++) {&#x2424; double[] w_3 = remove(model.sol(&quot;sol_4&quot;).getU(), 0);&#x2424; double[] w_2 = remove(model.sol(&quot;sol_3&quot;).getU(), 0);&#x2424; double[] w_1 = remove(model.sol(&quot;sol_2&quot;).getU(), 0);&#x2424; double[] w_0 = remove(model.sol(&quot;sol_1&quot;).getU(), 0);&#x2424; &#x2424; double p_3 = parametricFull[ii-1];&#x2424; double p_2 = parametricFull[ii-2];&#x2424; double p_1 = parametricFull[ii-3];&#x2424; double p_0 = parametricFull[ii-4];&#x2424; &#x2424; double[] W_3 = append(w_3, p_3);&#x2424; double[] W_2 = append(w_2, p_2);&#x2424; double[] W_1 = append(w_1, p_1);&#x2424; double[] W_0 = append(w_0, p_0);&#x2424; &#x2424; double As3 = 0;&#x2424; double As2 = 0;&#x2424; double As1 = 0;&#x2424; for (int i = 0; i &lt; W_3.length; i++) {&#x2424; As3 = As3+(W_3[i]-W_2[i])*(W_3[i]-W_2[i]);&#x2424; As2 = As2+(W_2[i]-W_1[i])*(W_2[i]-W_1[i]);&#x2424; As1 = As1+(W_1[i]-W_0[i])*(W_1[i]-W_0[i]);&#x2424; }&#x2424; double A3 = Math.sqrt(As3);&#x2424; double A2 = Math.sqrt(As2);&#x2424; double A1 = Math.sqrt(As1);&#x2424; double deltaA = A3;&#x2424; double yw0 = 0.0;&#x2424; double yw1 = yw0+A1;&#x2424; double yw2 = yw1+A2;&#x2424; double yw3 = yw2+A3;&#x2424; double yw4 = yw3+deltaA;&#x2424; double t0 = ((yw4-yw1)/(yw0-yw1))*((yw4-yw2)/(yw0-yw2))*((yw4-yw3)/(yw0-yw3));&#x2424; double t1 = ((yw4-yw0)/(yw1-yw0))*((yw4-yw2)/(yw1-yw2))*((yw4-yw3)/(yw1-yw3));&#x2424; double t2 = ((yw4-yw0)/(yw2-yw0))*((yw4-yw1)/(yw2-yw1))*((yw4-yw3)/(yw2-yw3));&#x2424; double t3 = ((yw4-yw0)/(yw3-yw0))*((yw4-yw1)/(yw3-yw1))*((yw4-yw2)/(yw3-yw2));&#x2424; double[] w_ini = new double[W_3.length];&#x2424; for (int i = 0; i &lt; W_3.length; i++) {&#x2424; w_ini[i] = t0*W_0[i]+t1*W_1[i]+t2*W_2[i]+t3*W_3[i];&#x2424; }&#x2424; //Initial values&#x2424; double[] W_ini = new double[w_ini.length-1];&#x2424; for (int i = 0; i &lt; w_ini.length-1; i++) {&#x2424; W_ini[i] = w_ini[i];&#x2424; }&#x2424; //Current value for the load&#x2424; double p_ini = w_ini[w_ini.length-1];&#x2424; &#x2424; parametricFull[ii] = p_ini;&#x2424; double[] W_ini_up = insert(W_ini, 1, 0);&#x2424; &#x2424; //Setting up the initial guess&#x2424; model.sol(&quot;sol_in&quot;).setU(W_ini_up);&#x2424; model.sol(&quot;sol_in&quot;).createSolution();&#x2424; &#x2424; //Load update&#x2424; model.param().set(p, p_ini);&#x2424; &#x2424; //Computing the current solution&#x2424; model.sol(&quot;sol_c&quot;).runAll();&#x2424; &#x2424; //Creating a parametric solution based on the current solution&#x2424; model.sol(&quot;sol_p&quot;).setU(ii+1, model.sol(&quot;sol_c&quot;).getU());&#x2424; &#x2424; //Updating the four last converged iterations&#x2424; model.sol(&quot;sol_2&quot;).copySolution(&quot;sol_1&quot;);&#x2424; model.sol(&quot;sol_3&quot;).copySolution(&quot;sol_2&quot;);&#x2424; model.sol(&quot;sol_4&quot;).copySolution(&quot;sol_3&quot;);&#x2424; model.sol(&quot;sol_c&quot;).copySolution(&quot;sol_4&quot;);&#x2424; &#x2424; //Debug log&#x2424; debugLog(&quot;Increment &quot;+toString(ii));&#x2424; debugLog(&quot;The load is &quot;+toString(p_ini));&#x2424; &#x2424; //Evaluating the stop condition&#x2424; if (isStopConditionGiven) {&#x2424; if (codeUtil.stopCondition(globalVariableMin, globalVariableMax)) {&#x2424; debugLog(&quot;Stop condition fulfilled for load &quot;+toString(p_ini));&#x2424; pfinal = ii;&#x2424; sizeOfSolution = ii+1;&#x2424; }&#x2424; }&#x2424;}&#x2424;&#x2424;double[] parametric = new double[sizeOfSolution];&#x2424;for (int ii = 0; ii &lt; sizeOfSolution; ii++) {&#x2424; parametric[ii] = parametricFull[ii];&#x2424;}&#x2424;model.sol(&quot;sol_p&quot;).setPVals(parametric);&#x2424;model.sol(&quot;sol_p&quot;).createSolution();&#x2424;if (contains(model.result().dataset().tags(), &quot;dset_c&quot;))&#x2424; model.result().dataset().remove(&quot;dset_c&quot;);&#x2424;long endTime = timeStamp();&#x2424;long totalTime = (endTime-startTime);&#x2424;String timeString = formattedTime(totalTime, &quot;hr:min:sec&quot;);&#x2424;debugLog(&quot;The solution time is &quot;+timeString);&#x2424;
//    In the Settings window for Method, locate the Inputs and Output section.
//    Find the Inputs subsection.
//    Click Add.
//     seven times.
//    In the table, enter the following settings:
//    In the Home toolbar, click Model Builder.
//     to switch to the main desktop.
//    Call the method <l>postBucklingWithArclength</l> in order to run it.
//    In the Home toolbar, click Model Builder.
//    In the Developer toolbar, click Method Call and choose postBucklingWithArclength.
//    In the Model Builder window, under Global Definitions, click PostBucklingWithArclength 1.
//    In the Settings window for Method Call, type Run Postbuckling Study with Arclength Method in the Label text field.
//    In the Tag text field, type postBucklingWithArclength.
//    Before running the <l>Run Postbuckling Study with Arclength Method</l>, assign correct inputs to the method.
//    Locate the Inputs section.
//    In the Dimensionless Load Parameter text field, type p.
//    In the Dimensionless Load Stepsize text field, type 0.0003.
//    In the Number of Iterations text field, type 200.
//    Select the Activate Stop Condition checkbox.
//    In the Global Variable in Stop Condition text field, type w_center.
//    In the Upper Limit of Global Variable text field, type 0.035.
//    Click Run.
//    Started method call postBucklingWithArclength

    model.param().set("p", "0");

    model.study().create("std_in");
    model.study("std_in").label("Placeholder Study for Initial Guess");
    model.study("std_in").create("stat_in", "Stationary");
    model.study("std_in").feature("stat_in").set("geometricNonlinearity", true);

    model.sol().create("sol_in");
    model.sol("sol_in").createAutoSequence("std_in");
    model.sol("sol_in").attach("std_in");

    model.study().create("std_c");
    model.study("std_c").label("Postbuckling Study");
    model.study("std_c").create("stat_c", "Stationary");
    model.study("std_c").feature("stat_c").set("geometricNonlinearity", true);

    model.sol().create("sol_c");
    model.sol("sol_c").createAutoSequence("std_c");
    model.sol("sol_c").attach("std_c");

    model.study().create("std_p");
    model.study("std_p").label("Placeholder Study for Parametric Solution");
    model.study("std_p").create("stat_p", "Stationary");

    model.sol().create("sol_p");
    model.sol("sol_p").createAutoSequence("std_p");
    model.sol("sol_p").attach("std_p");

    model.study("std_p").feature("stat_p").set("geometricNonlinearity", true);
    model.study("std_p").feature("stat_p").set("useparam", true);
    model.study("std_p").feature("stat_p").setIndex("pname", "p", 0);
    model.study("std_p").feature("stat_p").setIndex("plistarr", "1 2", 0);
    model.study("std_p").feature("stat_p").setIndex("punit", "N", 0);

    model.result().dataset().create("dset_c", "Solution");
    model.result().dataset("dset_c").set("solution", "sol_c");
    model.result().evaluationGroup().create("eg_c", "EvaluationGroup");
    model.result().evaluationGroup("eg_c").set("data", "dset_c");
    model.result().evaluationGroup("eg_c").create("gev1", "EvalGlobal");
    model.result().evaluationGroup("eg_c").feature("gev1").setIndex("expr", "w_center", 0);
    model.result().evaluationGroup("eg_c").set("includeparameters", false);
    model.result().dataset().remove("dset_c");

//    Finished method call postBucklingWithArclength
//    In the Home toolbar, click Model Builder.
//    In the Results toolbar, click 3D Plot Group.

    model.result().create("pg1", "PlotGroup3D");
    model.result("pg1").run();

//    In the Settings window for 3D Plot Group, type Stress in the Label text field.

    model.result("pg1").label("Stress");

//    Locate the Data section.
//    From the Dataset list, select Placeholder Study for Parametric Solution/Solution 3 (sol_p).

    model.result("pg1").set("data", "dset3");

//    Right-click Stress and choose Surface.

    model.result("pg1").create("surf1", "Surface");
    model.result("pg1").feature("surf1").set("evaluationsettings", "parent");

//    In the Settings window for Surface, locate the Expression section.
//    In the Expression text field, type shell.misesGp.

    model.result("pg1").feature("surf1").set("expr", "shell.misesGp");

//    Locate the Coloring and Style section.
//    From the Color table list, select Prism.

    return model;
  }

  public static Model run2(Model model) {

    model.result("pg1").feature("surf1").set("colortable", "Prism");

//    Right-click Surface 1 and choose Deformation.

    model.result("pg1").feature("surf1").create("def1", "Deform");
    model.result("pg1").run();

//    In the Settings window for Deformation, locate the Scale section.
//    Select the Scale factor checkbox.

    model.result("pg1").feature("surf1").feature("def1").set("scaleactive", true);

//    In the associated text field, type 1.

    model.result("pg1").feature("surf1").feature("def1").set("scale", 1);

//    In the Results toolbar, click Evaluation Group.

    model.result().evaluationGroup().create("eg1", "EvaluationGroup");

//    In the Settings window for Evaluation Group, type Evaluation Group: Force vs. Displacement in the Label text field.

    model.result().evaluationGroup("eg1").label("Evaluation Group: Force vs. Displacement");

//    Locate the Data section.
//    From the Dataset list, select Placeholder Study for Parametric Solution/Solution 3 (sol_p).

    model.result().evaluationGroup("eg1").set("data", "dset3");

//    Right-click Evaluation Group: Force vs. Displacement and choose Point Evaluation.

    model.result().evaluationGroup("eg1").create("pev1", "EvalPoint");

//    Select Point 4.

    model.result().evaluationGroup("eg1").feature("pev1").selection().set(4);

//    In the Settings window for Point Evaluation, locate the Expressions section.
//    In the table, enter the following settings:

    model.result().evaluationGroup("eg1").feature("pev1").setIndex("expr", "w_center", 0);
    model.result().evaluationGroup("eg1").feature("pev1").setIndex("expr", "P", 1);

//    In the Model Builder window, click Evaluation Group: Force vs. Displacement.
//    In the Settings window for Evaluation Group, click to expand the Format section.
//    From the Include parameters list, select Off.

    model.result().evaluationGroup("eg1").set("includeparameters", false);

//    In the Evaluation Group: Force vs. Displacement toolbar, click Evaluate.

    model.result().evaluationGroup("eg1").run();

//    In the Results toolbar, click 1D Plot Group.

    model.result().create("pg2", "PlotGroup1D");
    model.result("pg2").run();

//    In the Settings window for 1D Plot Group, type Force vs. Displacement in the Label text field.

    model.result("pg2").label("Force vs. Displacement");

//    Locate the Plot Settings section.
//    Select the x-axis label checkbox.

    model.result("pg2").set("xlabelactive", true);

//    In the associated text field, type Vertical displacement at shell center (m).

    model.result("pg2").set("xlabel", "Vertical displacement at shell center (m)");

//    Select the y-axis label checkbox.

    model.result("pg2").set("ylabelactive", true);

//    In the associated text field, type Force at shell center (N).

    model.result("pg2").set("ylabel", "Force at shell center (N)");

//    Locate the Legend section.
//    From the Position list, select Upper left.

    model.result("pg2").set("legendpos", "upperleft");

//    Right-click Force vs. Displacement and choose Table Graph.

    model.result("pg2").create("tblp1", "Table");
    model.result("pg2").feature("tblp1").set("markerpos", "datapoints");
    model.result("pg2").feature("tblp1").set("linewidth", "preference");

//    In the Settings window for Table Graph, locate the Data section.
//    From the Source list, select Evaluation group.

    model.result("pg2").feature("tblp1").set("source", "evaluationgroup");

//    Locate the Coloring and Style section.
//    Find the Line markers subsection.
//    From the Marker list, select Circle.

    model.result("pg2").feature("tblp1").set("linemarker", "circle");

//    From the Positioning list, select Interpolated.

    model.result("pg2").feature("tblp1").set("markerpos", "interp");

//    Click to expand the Legends section.
//    Select the Show legends checkbox.

    model.result("pg2").feature("tblp1").set("legend", true);

//    From the Legends list, select Manual.

    model.result("pg2").feature("tblp1").set("legendmethod", "manual");

//    In the table, enter the following settings:

    model.result("pg2").feature("tblp1").setIndex("legends", "Arc Length Approach", 0);

//    Right-click Table Graph 1 and choose Duplicate.

    model.result("pg2").feature().duplicate("tblp2", "tblp1");
    model.result("pg2").run();

//    In the Settings window for Table Graph, locate the Data section.
//    From the Source list, select Table.

    model.result("pg2").feature("tblp2").set("source", "table");

//    Locate the Coloring and Style section.
//    Find the Line markers subsection.
//    From the Marker list, select Diamond.

    model.result("pg2").feature("tblp2").set("linemarker", "diamond");

//    In the Number text field, type 12.

    model.result("pg2").feature("tblp2").set("markers", 12);

//    Locate the Legends section.
//    In the table, enter the following settings:

    model.result("pg2").feature("tblp2").setIndex("legends", "Global Equation Approach", 0);
    model.result("pg2").run();

//    In the Model Builder window, click Force vs. Displacement.
//    In the Force vs. Displacement toolbar, click Plot.

    model.result("pg2").run();

    model.title("Postbuckling Analysis Using an Incremental Arc Length Method");

    model
         .description("For slender structures, buckling is a catastrophic instability if the service load is above the critical limit. For such structures, it can be important to study the behavior of the structure beyond the critical buckling load, which is known as postbuckling analysis. Tracing the equilibrium path in postbuckling analysis is not easy as it leads to numerical difficulties such as limit points. Using an arc length method is a well-known strategy to trace equilibrium paths in such situations.\n\nIn this example, an incremental arc length method combined with cubic extrapolation is used for postbuckling analysis of a hinged cylindrical panel subjected to a point load at its center. The results of this example are compared with a similar example in which a global equation approach based on monotonically increasing displacement is used.");

    return model;
  }

  public static void main(String[] args) {
    Model model = run();
    run2(model);
  }

}
