4
$\begingroup$

I have an issue about a .wls script that I can't manage to solve for the moment. I execute it with wolframscript interpeter on command line :

A function PlotWithGrid has been coded this way :

PlotWithGrid[func_,initValueList_?ListQ, intervalList_?ListQ, axesLabels_ , PlotType_, minVec_, exportDir_?StringQ] := Module[
{ tb, f,cons, m, n, pl},
Switch[ PlotType,
1,
f[m_,n_, i_, j_]:=Apply[func,initValueList/.{initValueList[[i]]-> m, initValueList[[j]]-> n}];
tb = ParallelTable[
Piecewise[{{ContourPlot[f[x, y, i, j],{x,intervalList[[j,1]],intervalList[[j,2]]},{y,intervalList[[i,1]],intervalList[[i,2]]},AspectRatio->1,Frame->True,Axes->False, FrameLabel->{axesLabels[[j]], axesLabels[[i]]}],j<i}},Null],
{i, 2, Length@initValueList },
{j, 1, Length@initValueList - 1}
],
2,
tb = ParallelTable[
Piecewise[{{
ListPlot[func[[All, {j, i}]],Frame->True,FrameLabel->{axesLabels[[j]], axesLabels[[i]]},PlotStyle->{Red,PointSize[0.003]},PlotRange->{{intervalList[[j,1]],intervalList[[j,2]]}, {intervalList[[i,1]],intervalList[[i,2]]}},
Axes->False,ImageSize->Medium,
FrameTicks->{{Automatic, None}, {fticks, None}},
Epilog->{{PointSize[Medium],Blue,Point[func[[1,{j, i}]]],{PointSize[Medium],Black,Point[minVec[[{j, i}]]]}}},
AspectRatio->1
]
,j<i}},Null],
{i, 2, Length@initValueList },
{j, 1, Length@initValueList - 1}
],
3,
tb = Table[
Piecewise[{{
Show[
Confid2D[func[[3,All,{j, i}]],Blue],
Confid2D[func[[2,All,{j, i}]],Green],
Confid2D[func[[1,All,{j, i}]],Red],
PlotRange->{{intervalList[[j,1]],intervalList[[j,2]]}, {intervalList[[i,1]],intervalList[[i,2]]}},
FrameTicks->{{Automatic, None}, {fticks, None}},
Epilog->{PointSize[Medium],Black,Point[minVec[[{j, i}]]]},
FrameLabel->{axesLabels[[j]], axesLabels[[i]]},
AspectRatio->1,
ImageSize-> Medium
]
,j<i}},Null],
{i, 2, Length@initValueList },
{j, 1, Length@initValueList - 1}
],
4,
tb = Table[
Piecewise[{{
Show[
ListPlot[func[[3,All,{j, i}]], PlotStyle->{Blue,PointSize[Medium]}, Axes->False, Frame->True],
ListPlot[func[[2,All,{j, i}]], PlotStyle->{Green,PointSize[Medium]}, Axes->False, Frame->True],
ListPlot[func[[1,All,{j, i}]], PlotStyle->{Red,PointSize[Medium]}, Axes->False, Frame->True],
PlotRange->{{intervalList[[j,1]],intervalList[[j,2]]}, {intervalList[[i,1]],intervalList[[i,2]]}},
FrameTicks->{{Automatic, None}, {fticks, None}},
Epilog->{PointSize[Medium],Black,Point[minVec[[{j, i}]]]},
FrameLabel->{axesLabels[[j]], axesLabels[[i]]},
AspectRatio->1,
ImageSize-> Medium
]
,j<i}},Null],
{i, 2, Length@initValueList },
{j, 1, Length@initValueList - 1}
],
5,
cons[i_, j_]:=Quantile[PDF[SmoothKernelDistribution[func[[All,{ i, j}]]],func[[All, {i, j}]]],1-Erf[{1,2}/Sqrt[2]]//N];
f[ m_, n_, i_, j_]:=PDF[SmoothKernelDistribution[func[[All, {i, j}]]],{m, n}];

tb = Table[
Piecewise[{{
f[ m_, n_] = f[ m, n, j, i];
ContourPlot[f[m, n],{m, intervalList[[j, 1]],intervalList[[j,2]]},{n, intervalList[[i, 1]],intervalList[[i,2]]},Contours->cons[i, j],ContourShading->{White,Hue[0.5,0.9,0.9],Hue[0.55,0.8,0.8]},FrameLabel->{axesLabels[[j]], axesLabels[[i]]},PlotRange->{{intervalList[[j,1]],intervalList[[j,2]]}, {intervalList[[i,1]],intervalList[[i,2]]}},ClippingStyle->None, 
FrameTicks->{{Automatic, None}, {fticks, None}},Epilog->{PointSize[Medium],Black,Point[minVec[[{j, i}]]]}]
,j<i}},Null],
{i, 2, Length@initValueList },
{j, 1, Length@initValueList - 1}
]
];

pl=ResourceFunction["PlotGrid"][
tb,
PlotRange-> Max,
ImageSize-> Large,
Method->{"AllCustomTicks"->Automatic, "FixFrameTicks"->True}
];
Print@pl;
If[exportDir != "",
Export[exportDir,pl,"AllowRasterization"->False];
];

Select[tb//Flatten, #=!=Null&]

];
  1. If I use the function PlotWithGrid under this form :

    pl = PlotWithGrid[{data1s[[All, variableIdx]], data2s[[All, variableIdx]], data3s[[All, variableIdx]]}, {[CapitalOmega]m,[CapitalOmega]k,H0,[Phi]0,d[Phi]0,[Omega]BD}, varIntervals, axesLabels, 3,chainShort[[minimumChiPos]], "contour_1.png"];

The figure "contour_1.png" is empty and I don't understand why.

  1. I have into my code another type of plot like this :

    pl=ListPlot[chain[[All,8]],Joined->True,PlotRange->All,Frame->True,FrameLabel->{"N","!(*SuperscriptBox[([Chi]), (2)])"}]; Print@pl; pics = Join[pics, {pl}]; Export["chi2.png",pl,"AllowRasterization"->False];

In this case, the figure "chi2.png" is well generated (not empty).

  1. So I tried to use this Export function to my initial issue by doing :

    pl = PlotWithGrid[{data1s[[All, variableIdx]], data2s[[All, variableIdx]], data3s[[All, variableIdx]]}, {[CapitalOmega]m,[CapitalOmega]k,H0,[Phi]0,d[Phi]0,[Omega]BD}, varIntervals, axesLabels, 3,chainShort[[minimumChiPos]]]; pics = Join[pics, pl]; Export["contour_1.png",pl,"AllowRasterization"->False];

But with this snippet code, none figure "contour_1.png" is produced, nothing during the run.

If someone could indicate me a workaround to create correctly a none empty figure "contour_1.png" figure, this would be fine.

EDIT 1: Here a short example of .wls that don't work in direct command line : $ wolframswript short.wls

#!/usr/bin/env wolframscript
(* ::Package:: *)

(* SetOptions[EvaluationNotebook[],{AutoGeneratedPackage->Automatic,InitializationCellEvaluation->False,InitializationCellWarning->False,StyleDefinitions->Notebook[{Cell[StyleData[StyleDefinitions->"Default.nb"]],Cell[StyleData["Input"],InitializationCell->True]},Visible->False,StyleDefinitions->"PrivateStylesheetFormatting.nb"]}] *)


(* SetDirectory[NotebookDirectory[]]; *)
SetDirectory["."];
Needs["ComputationalGeometry`"];

fticks[min_,max_]:={#,Rotate[#,Pi/2],{ 0.01, 0}}&/@(N@FindDivisions[{min,max},6]);
PlotWithGrid[func_,initValueList_?ListQ, intervalList_?ListQ, axesLabels_ , PlotType_, minVec_, exportDir_?StringQ] := Module[
{ tb, f,cons, m, n, pl},
Switch[ PlotType,
2,
tb = ParallelTable[
Piecewise[{{
ListPlot[func[[All, {j, i}]],Frame->True,FrameLabel->{axesLabels[[j]], axesLabels[[i]]},PlotStyle->{Red,PointSize[0.003]},PlotRange->{{intervalList[[j,1]],intervalList[[j,2]]}, {intervalList[[i,1]],intervalList[[i,2]]}},
Axes->False,ImageSize->Medium,
FrameTicks->{{Automatic, None}, {fticks, None}},
Epilog->{{PointSize[Medium],Blue,Point[func[[1,{j, i}]]],{PointSize[Medium],Black,Point[minVec[[{j, i}]]]}}},
AspectRatio->1
]
,j<i}},Null],
{i, 2, Length@initValueList },
{j, 1, Length@initValueList - 1}
],
3,
tb = Table[
Piecewise[{{
Show[
Confid2D[func[[3,All,{j, i}]],Blue],
Confid2D[func[[2,All,{j, i}]],Green],
Confid2D[func[[1,All,{j, i}]],Red],
PlotRange->{{intervalList[[j,1]],intervalList[[j,2]]}, {intervalList[[i,1]],intervalList[[i,2]]}},
FrameTicks->{{Automatic, None}, {fticks, None}},
Epilog->{PointSize[Medium],Black,Point[minVec[[{j, i}]]]},
FrameLabel->{axesLabels[[j]], axesLabels[[i]]},
AspectRatio->1,
ImageSize-> Medium
]
,j<i}},Null],
{i, 2, Length@initValueList },
{j, 1, Length@initValueList - 1}
]
];

pl=ResourceFunction["PlotGrid"][
tb,
PlotRange-> Max,
ImageSize-> Large,
Method->{"AllCustomTicks"->Automatic, "FixFrameTicks"->True}
];
If[exportDir != "",
Export[exportDir,pl,"AllowRasterization"->False];
];

Select[tb//Flatten, #=!=Null&]

];

Confid2D[data_?ListQ,color_]:=Module[{cv,cv1,cv2},
cv=ConvexHull[data,AllPoints->True];
cv1=data[[#]]&/@cv;
If[Length@cv1 > 1, 
cv2=Join[cv1,{cv1[[1]],cv1[[2]]}],
cv2=Join[cv1,{cv1[[1]]}]];
Show[Graphics[{color,FilledCurve[BSplineCurve[cv2,SplineDegree->2,SplineKnots->"Unclamped"]],Antialiasing->True}],Frame->True]
];

dH[\[Rho]m_,\[Phi]_,u_,\[Omega]BD_, \[CapitalOmega]k_,du_, z_]:= -(((1+z) (16 N[\[Pi]] \[Rho]m-6 (1+z)^2 \[CapitalOmega]k \[Phi]) ((1+z) \[Omega]BD u^3-2 \[Omega]BD u ((1+z) du[\[Rho]m,\[Phi],u,\[Omega]BD, \[CapitalOmega]k,z]+u) \[Phi]-6 du[\[Rho]m,\[Phi],u,\[Omega]BD, \[CapitalOmega]k,z] \[Phi]^2)+(1/(1+z))6 \[Phi] (-8 N[\[Pi]] \[Rho]m+(1+z)^2 \[CapitalOmega]k ((1+z) u+2 \[Phi])) (6 \[Phi]^2-(1+z) u ((1+z) \[Omega]BD u+6 \[Phi])))/(2 Sqrt[(-16 N[\[Pi]] \[Rho]m+6 (1+z)^2 \[CapitalOmega]k \[Phi])/(6 (1+z) u+((1+z)^2 \[Omega]BD u^2)/\[Phi]-6 \[Phi])] ((1+z)^2 \[Omega]BD u^2+6 (1+z) u \[Phi]-6 \[Phi]^2)^2));
d\[Rho]m[\[Rho]m_,\[Phi]_,u_,z_]:=3/(1+z) \[Rho]m;
d\[Phi][\[Rho]m_,\[Phi]_,u_,z_]:= u;
du[\[Rho]m_,\[Phi]_,u_,\[Omega]BD_, \[CapitalOmega]k_, z_]:=-((24 N[\[Pi]] \[Rho]m \[Phi]^3+(1+z) u \[Phi]^2 (8 \[Pi] (-3+\[Omega]BD) \[Rho]m-3 (1+z)^2 (3+2 \[Omega]BD) \[CapitalOmega]k \[Phi])+3 (1+z)^2 u^2 \[Phi] (-4 N[\[Pi]] \[Omega]BD \[Rho]m+(1+z)^2 (3+2 \[Omega]BD) \[CapitalOmega]k \[Phi])+(1+z)^3 \[Omega]BD u^3 (-4 \[Pi] (1+\[Omega]BD) \[Rho]m+(1+z)^2 (3+2 \[Omega]BD) \[CapitalOmega]k \[Phi]))/((1+z)^2 (3+2 \[Omega]BD) \[Phi]^2 (-8 N[\[Pi]] \[Rho]m+3 (1+z)^2 \[CapitalOmega]k \[Phi])));

RK4Method[dH_,d\[Phi]_,d\[Rho]m_,du_,\[CapitalOmega]m_,\[CapitalOmega]k_,H0_,\[Phi]0_,d\[Phi]0_,\[Omega]BD_,zLine0_]:=Module[
{ h, Htable, \[Rho]mtable, \[Phi]table, utable, Hk1, Hk2, Hk3,Hk4, \[Rho]mk1, \[Rho]mk2, \[Rho]mk3,\[Rho]mk4, \[Phi]k1, \[Phi]k2,\[Phi]k3,\[Phi]k4, uk1,uk2, uk3,uk4, containsIndeterminate, containsComplex, Hval,zLine =Join[{0},zLine0]},
containsIndeterminate = False;
containsComplex = False;

Htable=Table[Infinity,{i,1,Length@zLine}];
Htable[[1]]=H0;
\[Rho]mtable=3 H0 H0 \[Phi]0 \[CapitalOmega]m/(8N[\[Pi]]);
utable=d\[Phi]0;
\[Phi]table=\[Phi]0;
Do[
h=(zLine[[i]]-zLine[[i - 1]]);
Hk1=dH[\[Rho]mtable,\[Phi]table,utable,\[Omega]BD, \[CapitalOmega]k,du,zLine[[i - 1]]];
\[Phi]k1=d\[Phi][\[Rho]mtable,\[Phi]table,utable,zLine[[i - 1]]];
\[Rho]mk1=d\[Rho]m[\[Rho]mtable,\[Phi]table,utable,zLine[[i - 1]]];
uk1=du[\[Rho]mtable,\[Phi]table,utable,\[Omega]BD, \[CapitalOmega]k,zLine[[i - 1]]];

Hk2=dH[\[Rho]mtable+h/2*\[Rho]mk1, \[Phi]table+h/2* \[Phi]k1,utable + h/2* uk1, \[Omega]BD, \[CapitalOmega]k,du, zLine[[i - 1]]+h/2];
\[Phi]k2=d\[Phi][\[Rho]mtable+h/2*\[Rho]mk1, \[Phi]table+h/2* \[Phi]k1,utable + h/2* uk1,zLine[[i - 1]]+h/2];
\[Rho]mk2=d\[Rho]m[\[Rho]mtable+h/2*\[Rho]mk1, \[Phi]table+h/2* \[Phi]k1,utable + h/2* uk1,zLine[[i - 1]]+h/2];
uk2=du[\[Rho]mtable+h/2*\[Rho]mk1, \[Phi]table+h/2* \[Phi]k1,utable+ h/2* uk1,\[Omega]BD, \[CapitalOmega]k,zLine[[i - 1]]+h/2];

Hk3=dH[\[Rho]mtable+h/2*\[Rho]mk2,\[Phi]table+h/2* \[Phi]k2,utable+h/2* uk2,\[Omega]BD, \[CapitalOmega]k,du,zLine[[i - 1]]+h/2];
\[Phi]k3=d\[Phi][\[Rho]mtable+h/2*\[Rho]mk2,\[Phi]table+h/2* \[Phi]k2,utable+h/2* uk2,zLine[[i - 1]]+h/2];
\[Rho]mk3=d\[Rho]m[\[Rho]mtable+h/2*\[Rho]mk2,\[Phi]table+h/2* \[Phi]k2,utable+h/2* uk2,zLine[[i - 1]]+h/2];
uk3=du[\[Rho]mtable+h/2*\[Rho]mk2,\[Phi]table+h/2* \[Phi]k2,utable+h/2* uk2,\[Omega]BD, \[CapitalOmega]k,zLine[[i - 1]]+h/2];

Hk4=dH[\[Rho]mtable+h*\[Rho]mk3,\[Phi]table+h* \[Phi]k3,utable+h* uk3,\[Omega]BD, \[CapitalOmega]k,du,zLine[[i ]]];
\[Phi]k4=d\[Phi][\[Rho]mtable+h*\[Rho]mk3,\[Phi]table+h* \[Phi]k3,utable+h* uk3,zLine[[i ]]];
\[Rho]mk4=d\[Rho]m[\[Rho]mtable+h*\[Rho]mk3,\[Phi]table+h* \[Phi]k3,utable+h* uk3,zLine[[i ]]];
uk4=du[\[Rho]mtable+h*\[Rho]mk3,\[Phi]table+h* \[Phi]k3,utable+h* uk3,\[Omega]BD, \[CapitalOmega]k,zLine[[i ]]];

Hval=Htable[[i-1]]+h (Hk1+2 Hk2+2 Hk3+Hk4)/6;

If[
Hval===Indeterminate,
containsIndeterminate = True;
Break[]
];

If[
!Element[Hval, Reals],
containsComplex = True;
Break[]
];
Htable[[i]] = Hval;
\[Rho]mtable=\[Rho]mtable+h (\[Rho]mk1+2 \[Rho]mk2+2 \[Rho]mk3+\[Rho]mk4)/6;
\[Phi]table=\[Phi]table+h (\[Phi]k1+2 \[Phi]k2+2 \[Phi]k3+\[Phi]k4)/6;
utable=utable+h (uk1+2 uk2+2 uk3+uk4)/6,

{i,2,Length@zLine }
];
{containsIndeterminate, containsComplex, Htable[[2;;]]}];

applyMCMC[variables_?ListQ, chi2old_, priors_?ListQ]:=Module[{\[CapitalOmega]m,\[CapitalOmega]k,H0,\[Phi]0,d\[Phi]0,\[Omega]BD, chi2prop, \[CapitalOmega]mNew,\[CapitalOmega]kNew,\[CapitalOmega]deNew,H0New,\[Phi]0New,d\[Phi]0New,\[Omega]BDNew, containsIndeterminate, containsComplex, vals, jumpprob, newVariables},
(* Propose a new point based on the jump function with different standard dev*)

{\[CapitalOmega]mNew,\[CapitalOmega]kNew,H0New,\[Phi]0New,d\[Phi]0New,\[Omega]BDNew} = variables +( RandomVariate[NormalDistribution[0,#]]&/@{0.01,0.001, 0.1,0.1, 0.001, 10 });

(*{\[CapitalOmega]mNew,\[CapitalOmega]kNew,H0New,\[Phi]0New,d\[Phi]0New,\[Omega]BDNew} = RandomVariate[MultinormalDistribution[variables,DiagonalMatrix[{0.01,0.01, 0.1,0.01, 0.001, 1.5 }]^2]]*);

\[CapitalOmega]deNew=1 - \[CapitalOmega]mNew - \[CapitalOmega]kNew;

{containsIndeterminate, containsComplex, vals} = Quiet@RK4Method[dH,d\[Phi],d\[Rho]m,du,\[CapitalOmega]mNew,\[CapitalOmega]kNew,H0New,\[Phi]0New,d\[Phi]0New,\[Omega]BDNew,zLine];

newVariables = {};

If[!containsIndeterminate ,

If[!containsComplex ,
(* This is the \[Chi]^2 at the proposed point *)
chi2prop=chi2[vals];
(* This is the probability to jump to the new point=Exp[-\[Chi]^2(prop)]/Exp[-\[Chi]^2(current)]*)
jumpprob=Exp[-(chi2prop-chi2old)/2];
(*Print@jumpprob;*)
(* Decide whether to keep the new point taking into account also priors. If so, store it in the chain, otherwise start over. *)
If[RandomReal[]<=Min[{1,jumpprob}]&& And@@ IntervalMemberQ[Interval@#&/@priors,{\[CapitalOmega]mNew,\[CapitalOmega]kNew,H0New,\[Phi]0New,d\[Phi]0New,\[Omega]BDNew}],
newVariables = {\[CapitalOmega]mNew,\[CapitalOmega]kNew,\[CapitalOmega]deNew,H0New,\[Phi]0New,d\[Phi]0New,\[Omega]BDNew,chi2prop};
];
];
];
newVariables
];

dataH=Import["H_All.txt","Table"];
dataH = DeleteCases[dataH, x_?(Length[#]==0&), 1];
ndata=Length[dataH];
zLine = dataH[[All, 1]];

(* Choose the total number of steps *)
totalsteps=30000;
(* Choose the initial point for BD osmology as {\[CapitalOmega]m,\[CapitalOmega]de,\[CapitalOmega]k=1.0-\[CapitalOmega]m-\[CapitalOmega]de,H0,\[Phi]0,d\[Phi]0,\[Omega]BD}*)
\[CapitalOmega]m=0.275;
\[CapitalOmega]k=-0.005;
\[CapitalOmega]de=1-\[CapitalOmega]m - \[CapitalOmega]k;
H0=71.0;
\[Phi]0=2.5;
d\[Phi]0=0.005;
\[Omega]BD=77.75;
burnin=500;
priors={{0.24,0.4},{-0.01,0.01}, {64,76}, {0,4}, {0., 0.01}, {-350,350}};
variables = {\[CapitalOmega]m,\[CapitalOmega]k,H0,\[Phi]0,d\[Phi]0,\[Omega]BD};
variableIdx = {1,2,4,5, 6, 7};
axesLabels = {"\[CapitalOmega]m","\[CapitalOmega]k","H0","\[Phi]0","d\[Phi]0","\[Omega]BD"};
chi2[h_?ListQ]:=Total[((dataH[[All,2]]-h)/dataH[[All,3]])^2];

(* This is the current point solution*)
{containsIndeterminate, containsComplex, vals}  = RK4Method[dH,d\[Phi],d\[Rho]m,du,\[CapitalOmega]m,\[CapitalOmega]k,H0,\[Phi]0,d\[Phi]0,\[Omega]BD,zLine];
(* This is the \[Chi]^2 at the current point*)
chi2old= chi2[vals];

(* Initialize the chain. This is were the points are stored *)
chain={};
chlen = 0;

(* Start the MCMC *)
(* Initialize the MCMC *)

step=0;
While[
step< totalsteps,
newVariables=Quiet@applyMCMC[variables, chi2old, priors];
If[Length@newVariables>1,
AppendTo[chain,newVariables];
variables = newVariables[[variableIdx]];
chi2old = newVariables[[Length@newVariables]];
chlen=Dimensions[chain][[1]];
];
step++
]

Export["chain.m",chain];
(* The values of \[CapitalDelta]\[Chi]^2 for 1,2,3\[Sigma] *)
dx2=N@Table[2 InverseGammaRegularized[1/2,1-Erf[x\[Sigma]/(\[Sqrt]2)]],{x\[Sigma],1,3}];
chi2List = chain[[All, Dimensions[chain][[2]]]];
Print[StringForm["Estimated Initial \!\(\*SuperscriptBox[\(\[Chi]\), \(2\)]\) is ``, \nEstimated Minimum \!\(\*SuperscriptBox[\(\[Chi]\), \(2\)]\) is ``",chi2List[[1]],Min@chi2List]];
(*sortedIdx = Reverse@Ordering[chi2List];
chain = chain[[sortedIdx, All]];*)
chi2min = Min@chi2List;
minimumChiPos = Flatten[ Position[chi2List, chi2min], 1]//First;
{\[CapitalOmega]m,\[CapitalOmega]k,\[CapitalOmega]de,H0,\[Phi]0,d\[Phi]0,\[Omega]BD,chi2current} = chain[[minimumChiPos, All]];
Print[StringForm["Parameters converged to the following values \n\[CapitalOmega]m -> ``,\n \[CapitalOmega]k -> ``,\n \[CapitalOmega]de -> ``,\n H0 -> ``,\n \[Phi]0 -> ``,\n d\[Phi]0 -> ``,\n \[Omega]BD -> ``",\[CapitalOmega]m,\[CapitalOmega]k,\[CapitalOmega]de,H0,\[Phi]0,d\[Phi]0,\[Omega]BD]];

chainShort = chain[[All,variableIdx ]];
varIntervals = MinMax/@Transpose@chainShort;

pl = PlotWithGrid[chainShort, {\[CapitalOmega]m,\[CapitalOmega]k,H0,\[Phi]0,d\[Phi]0,\[Omega]BD}, varIntervals, axesLabels, 2,chainShort[[minimumChiPos]], "scatter.png"];
(* Select the points in the 1,2,3 sigma region *)
data1s=Select[chain,#[[8]]<=(chi2min+dx2[[1]])&];
data2s=Select[chain,#[[8]]<=(chi2min+dx2[[2]])&];
data3s=Select[chain,#[[8]]<=(chi2min+dx2[[3]])&];
varIntervals = MinMax/@Transpose@data3s[[All, variableIdx]];

pl = PlotWithGrid[{data1s[[All, variableIdx]], data2s[[All, variableIdx]], data3s[[All, variableIdx]]}, {\[CapitalOmega]m,\[CapitalOmega]k,H0,\[Phi]0,d\[Phi]0,\[Omega]BD}, varIntervals, axesLabels, 3,chainShort[[minimumChiPos]], "contour_1.png"];

**Results : scatter.png  and contour_1.png are empty.**

The file "H_All.txt" is below :

0.0700 69.00 19.6
0.0900 69.00 12.0
0.1200 68.60 26.2
0.1700 83.00 8.00
0.1790 75.00 4.00
0.1990 75.00 5.00
0.2000 72.90 29.6
0.2400 79.69 2.65
0.2700 77.00 14.0
0.2800 88.80 36.6
0.3500 84.40 7.00
0.3520 83.00 14.0
0.3800 83.00 13.5
0.4000 95.00 17.0
0.4004 77.00 10.2
0.4247 87.10 11.2
0.4300 86.45 3.68
0.4400 82.60 7.80
0.4449 92.80 12.9
0.4783 80.90 9.00
0.4800 97.00 62.0
0.5700 92.40 4.50
0.5930 104.0 13.0
0.6000 87.90 6.10
0.6800 92.00 8.00
0.7300 97.30 7.00
0.7810 105.0 12.0
0.8750 125.0 17.0
0.8800 90.00 40.0
0.9000 117.0 23.0
1.0370 154.0 20.0
1.3000 168.0 17.0
1.3630 160.0 33.6
1.4300 177.0 18.0
1.5300 140.0 14.0
1.7500 202.0 40.0
1.9650 186.5 50.4
2.3000 224.0 8.00
2.3400 222.0 7.00
2.3600 226.0 8.00
    

and the error log I get during a direct launching by command line :

$ wolframscript BD_MCMC_2_short.wls

StringForm[Estimated Initial DisplayForm[SuperscriptBox[χ, 2]] is ``,
Estimated Minimum DisplayForm[SuperscriptBox[χ, 2]] is ``, 220.3044990025602, 21.58177375066642]
StringForm[Parameters converged to the following values
Ωm -> ``,
 Ωk -> ``,
 Ωde -> ``,
 H0 -> ``,
 ϕ0 -> ``,
 dϕ0 -> ``,
 ωBD -> ``, 0.24041269047499136, 0.0016398283879296882, 0.757947481137079, 64.06857800452705, 1.579419511447762, 0.00011446222576435193, 334.91426614586453]

FrontEndObject::notavail: A front end is not available; certain operations require a front end.

FrontEndObject::notavail: A front end is not available; certain operations require a front end.

FrontEndObject::notavail: A front end is not available; certain operations require a front end.

General::stop: Further output of FrontEndObject::notavail will be suppressed during this calculation.

ReplaceAll::reps: {$Failed} is neither a list of replacement rules nor a valid dispatch table, and so cannot be used for replacing.

Symbol::argx: Symbol called with 0 arguments; 1 argument is expected.

Symbol::argx: Symbol called with 0 arguments; 1 argument is expected.

Function::slotn: Slot number 2 in {Reverse[MinMax[#1]], MinMax[#2]} &  cannot be filled from ({Reverse[MinMax[#1]], MinMax[#2]} & )[Missing[String[], Symbol[]]].

Function::slotn: Slot number 2 in {Reverse[MinMax[#1]], MinMax[#2]} &  cannot be filled from ({Reverse[MinMax[#1]], MinMax[#2]} & )[String[]].

Function::slotn: Slot number 2 in {Reverse[MinMax[#1]], MinMax[#2]} &  cannot be filled from ({Reverse[MinMax[#1]], MinMax[#2]} & )[Symbol[]].

General::stop: Further output of Function::slotn will be suppressed during this calculation.

Outer::heads: Heads List and Missing at positions 3 and 2 are expected to be the same.

...                                                                                                                                                                                                                Total[Transpose[Mean[MapThread[Switch[(Abs @* Subtract)[1], Inherited, (Abs @* Subtract)[2], Min, (Abs @* Subtract)[#3, 1], Max, (Abs @* Subtract)[#3, 2], (Abs @* Subtract)[], (Abs @* Subtract)[1]] & , {{{{Abs[Subtract[Max, Max, Max, Max, Max]], Abs[Subtract[Max, Max, Max, Max, Max]]}, {Abs[Subtract[Max, Max, Max, Max, Max]], Abs[Subtract[Max, Max, Max, Max, Max]]}}, {{Abs[Subtract[Max, Max, Max, Max, Max]], Abs[Subtract[Max, Max, Max, Max, Max]]}, {Abs[Subtract[Max, Max, Max, Max, Max]], Abs[Subtract[Max, Max, Max, Max, Max]]}}, {{Abs[Subtract[Max, Max, Max, Max, Max]], Abs[Subtract[Max, Max, Max, Max, Max]]}, {Abs[Subtract[Max, Max, Max, Max, Max]], Abs[Subtract[Max, Max, Max, Max, Max]]}}, {{Abs[Subtract[Max, Max, Max, Max, Max]], Abs[Subtract[Max, Max, Max, Max, Max]]}, {Abs[Subtract[Max, Max, Max, Max, Max]], Abs[Subtract[Max, Max, Max, Max, Max]]}}, {{Abs[Subtract[Max, Max, Max, Max, Max]], Abs[Subtract[Max, Max, Max, Max, Max]]}, {Abs[Subtract[Max, Max, Max, Max, Max]], Abs[Subtract[Max, Max, Max, Max, Max]]}}}, Missing[NotAvailable, PlotRange], {{{Abs[Subtract[#2]], Abs[Subtract[#1]]}[{Abs[String[] - Symbol[]], Abs[String[] - Symbol[]]}, 0], {Abs[Subtract[#2]], Abs[Subtract[#1]]}[{Abs[String[] - Symbol[]], Abs[String[] - Symbol[]]}, 0]}, {{Abs[Subtract[#2]], Abs[Subtract[#1]]}[{0, 0}, 0], {Abs[Subtract[#2]], Abs[Subtract[#1]]}[{0, 0}, 0]}}}, 4]]]]   Mean[1] + Mean[2] + Mean[3]  Mean[1] + Mean[2] + Mean[3]  Mean[1] + Mean[2] + Mean[3]              Total[Transpose[Mean[Missing[NotAvailable, ImageSize]]]]   Mean[1] + Mean[2] + Mean[3]  Mean[1] + Mean[2] + Mean[3]  Mean[1] + Mean[2] + Mean[3]
    
    Transpose::argt: Transpose called with 0 arguments; 1 or 2 arguments are expected.
    
    General::stop: Further output of Transpose::argt will be suppressed during this calculation.
    
    MapThread::list: List expected at position 2 in MapThread[If[#2 === Null, ##1] & , Missing[PartAbsent, 2], 2].
    
    Divide::argr: Divide called with 1 argument; 2 arguments are expected.
    
    ^C
    Interrupt> Your options are:
        abort (or a) to abort current calculation
        continue (or c) to continue
        exit (or quit) to exit Mathematica
        inspect (or i) to enter an interactive dialog

EDIT 2: Here a suggestion from Wolfram support :

" The appropriate syntax would be

UsingFrontEnd[pl=ResourceFunction["PlotGrid"][...]]

However, please note that ResourceFunction items are not part of the core Wolfram Language, and are instead created by community members. As such, we cannot provide assistance or guidance with any specific ResourceFunction, and any issues regarding PlotGrid not working with wolframscript or when run in certain ways should be directed to the function's authors via the feedback form "

So, i tried to do :

UsingFrontEnd[pl=ResourceFunction["PlotGrid"][
tb,
PlotRange-> Max,
ImageSize-> Large,
Method->{"AllCustomTicks"->Automatic, "FixFrameTicks"->True}
]];

But get still errors during execution of Wolframscript .wls. Any help is welcome.

$\endgroup$
6
  • $\begingroup$ Can you clarify how you are running your code exactly? You mention wolframscript, but you don't show the command that you are using, and you appear to be using data that are not given in your question. It would be a lot easier to diagnose your issues if you could tell us exactly what you're evaluating $\endgroup$
    – Lukas Lang
    May 15, 2022 at 20:25
  • $\begingroup$ @LukasLang . I am running the script "main_script.wls" like this : $ wolframscript main_script.wls . What is weird is that figures are generated under the form of Notebook and with previous direct command line cited, the figure are empty. I will add some additional informations about the code. However, I don't understand why it works with Notebook and not with terminal launching command. $\endgroup$
    – youpilat13
    May 15, 2022 at 21:21
  • $\begingroup$ Please include a fully working example in your question (you can make sure you can run the code on your computer by copying into a new file and trying to run that). This includes some example-data for data1s etc., otherwise we have to guess how their structure looks like $\endgroup$
    – Lukas Lang
    May 16, 2022 at 10:16
  • $\begingroup$ Thanks for providing the full code. Unfortunately, it appears that there is an issue where some of the functionality required by PlotGrid is not available from wolframscript. After trying around for a bit, I don't think this is something that can easily be fixed $\endgroup$
    – Lukas Lang
    May 17, 2022 at 9:37
  • $\begingroup$ @LukasLang I have added an EDIT 2 possible workaround but I get still errors during execution of .wls script. $\endgroup$
    – youpilat13
    May 17, 2022 at 22:16

1 Answer 1

5
+50
$\begingroup$

Ok, finally managed to get it to run from within wolframscript...

I have created an updated version of PlotGrid that should work properly in non-standard environments such as wolframscript. It is not yet published, but you can get it as ResourceFunction@CloudObject["https://www.wolframcloud.com/obj/langl/DeployedResources/Function/PlotGrid/ResourceObject.wl"] in the meantime.

With this, your script produces plots as expected:

(* rest of script ... *)

pl=ResourceFunction[CloudObject["https://www.wolframcloud.com/obj/langl/DeployedResources/Function/PlotGrid/ResourceObject.wl"]][
   tb,
   PlotRange-> Max,
   ImageSize-> Large,
   Method->{"AllCustomTicks"->Automatic, "FixFrameTicks"->True}
];

(* rest of script ... *)

The resulting contour_1.png and scatter.png look like this for me:

enter image description here

enter image description here

What went wrong

For those who are curious: There are a few things that went wrong in the original version of PlotGrid

  • UsingFrontEnd was required since PlotGrid uses CallFrontEnd internally to obtain the dimensions of the plots (see also this answer by @CarlWoll). I have now included UsingFrontEnd, so that it automatically starts a frond-end if required.
  • The main issues all come from the fact that the front-end launched by UsingFrontEnd in turn launches another kernel that is used to evaluate GridLines, Ticks and FrameTicks. This has two consequences:
    • Using custom ticks functions doesn't work, since the other kernel doesn't know about them. This is fixed rather easily by transporting the relevant definitions to the other kernel using Language`ExtendedFullDefinition inside a FrameTicks specification.

    • Any tracking injected into these options can't communicate the results back to the main kernel. To complicate things further, the other kernel is launched in a sandboxed mode, heavily restricting the communication channels. The best solution I have found is a trick based on ExportPacket:

      FrontEndExecute@FrontEnd`ExportPacket[
        Cell@BoxData@DynamicBox@ToBoxes[expr],
        "PlainText"
      ]
      

      This evaluates expr on the other kernel and returns the result as a string, which we can then interpret using ToExpression. Starting in version 13.0, this communication is no longer required, since we can get all the necessary information from the newly-fixed AbsoluteOptions. The definitions however still need to be transported over.

In total, the strategy is thus as follows:

  1. "Rasterize" a Graphics expression with FrameTicks->{{None,None},{None,(Language`ExtendedFullDefinition[]=def)&}. This transports the relevant definitions to the other kernel. Here, "rasterize" is effectively equivalent to Rasterize[...,"Regions"].
  2. "Rasterize"/Measure the actual Graphics expressions. Here, we save some data about the plot range in the DownValues of SowData on the helper kernel. Originally, this was using Sow, but since this stuff is now evaluated on the other kernel, it didn't make its way back to the main kernel.
  3. Use the ExportPacket trick to get the DownValues of SowData to the main kernel, and evaluate them, effectively replaying the Sow calls.
$\endgroup$
5
  • $\begingroup$ Thanks a lot ! could you send me please a link of the script, I am not confortable with the cloud of wolfram. Thanks in advance, best regards $\endgroup$
    – youpilat13
    May 23, 2022 at 22:25
  • $\begingroup$ You can download the definition notebook from here: dropbox.com/s/eu1k4br8nhj7jpo/PlotGrid-4-3-0.nb?dl=0. You can open the notebook, and press "Deploy > Locally to this computer" from the orange toolbar. Afterwards, it will take some minutes to process, and you should get a banner notification with something like "Deployed to your computer". You can then click the small button next to the "PlotGrid" icon in the toolbar (tooltip: "Click to copy to clipboard"). You should then be able to paste a reference to the resource function, something like ResourceFunction[LocalObject[...]]. $\endgroup$
    – Lukas Lang
    May 24, 2022 at 7:30
  • $\begingroup$ excuse me but could send me please directly the entire Wolfram script source .wls, I am not good to convert Mathematica Notebook to wolfram script. I would be grateful if you could it. Regards $\endgroup$
    – youpilat13
    May 24, 2022 at 9:04
  • $\begingroup$ You don't really have to convert the notebook to a .wls script: You just have to create the resource function once using the procedure I outline above, and then you can use the ResourceFunction[LocalObject[...]] expression directly in your .wls script without any issues. That being said, I tried to create a .wls version of the function here: dropbox.com/s/bp3hv9jlg50dx2u/PlotGrid-4-3-0.wls?dl=0. You use it as Get["path/to/file.wls"], and then the function is simply used like PlotGrid[tb,...]. $\endgroup$
    – Lukas Lang
    May 24, 2022 at 9:28
  • $\begingroup$ Thank you so much ! I am going to test it $\endgroup$
    – youpilat13
    May 24, 2022 at 10:57

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.