Skip to main content

plotting - Exporting the solutions of two ODEs sampled over same set of domain values



I have two functions h[r] and c[r] defined by two differentials equation.


operh[h_, c_] := h D[D[h, r, r] + 1/r D[h, r] - 1/lc^2 h, r] + m D[c, r]
operc[h_, c_] := 
  cmt h + ξ ((ϕ - H) (1 - ϕ) )/Sqrt[R^2 - r^2] - 
    1/r D[r h (1 + e^2/1680  h^4  (D[c, r])^2 ) D[c, r], r]

I numerically solve them with NDSolve (with solution called "s" here under). Then I plot them with no problem:



Plot[h[r] /. s, {r, rleft, rright}, PlotRange -> All, AxesOrigin -> {0, 0}]
Plot[c[r] /. s, {r, rleft, rright}, PlotRange -> All, AxesOrigin -> {0, 0}]

With same rright and rleft.


graph for h graph for c


What I want is a .txt or .csv file with 3 columns, containing r, h[r] and c[r]. From what I have read on this site I tried:


ploth = 
  Plot[h[r] /. s, {r, rleft, rright}, PlotRange -> All, AxesOrigin -> {0, 0}]
pts = Cases[ploth, Line[pts_] :> pts, Infinity][[1]]
Export[SystemDialogInput["FileSave", "filename.csv"], pts]


It gave me a list of 157 pairs {r, h[r]}. Then I did the same for c:


plotc = 
  Plot[c[r] /. s, {r, rleft, rright}, PlotRange -> All, AxesOrigin -> {0, 0}]
pts = Cases[plotc, Line[pts_] :> pts, Infinity][[1]]
Export[SystemDialogInput["FileSave", "filename.csv"], pts]

And it gave me a list of 226 pairs {r, c[r]}. This is not the same number of data than for h(r), so I can't merge them.


I tried by plotting the two curves on a single graph but it gave me only one of the two data sets.


I also tried to use Table but it gave me an empty list:



tableh = Table[{r, ploth}, {r, rleft, rright}]


{}



Do you know how to extract h and c data for common r values?


I'm on Mathematica 7.0



Answer



An approach to this issue I like is to rely on the fact that NDSolve found the solution on the same grid points for both functions in the first place:


a simple example:



 sol = First@NDSolve[ {f''[x] + g'[x] + 1 == 0 , g[x] == x^2, f[0] == 0, 
   f'[0] == 1}, {f, g}, {x, 0, 1}];
 Plot[{f[x] /. sol, (g[x] /. sol)}, {x, 0, 1}]

enter image description here


here we do not actually use the interpolation but extract the grid values:


 gridpoints = Flatten[(g /. sol)["Grid"]];
 computeddata = 
    Transpose[{gridpoints, 
              (f /. sol)["ValuesOnGrid"],

              (g /. sol)["ValuesOnGrid"]}];
 ListPlot[{ computeddata[[All, {1, 2}]], computeddata[[All, {1, 3}]]}, 
  Joined -> True, PlotMarkers -> Automatic]

enter image description here


here computeddata is a nice array you can simply export with each row {x,f[x],g[x]}


An advantage to this is the point spacing is naturally refined as was needed by NDSolve


Comments

Post a Comment

Popular posts from this blog

plotting - Plot 4D data with color as 4th dimension

I have a list of 4D data (x position, y position, amplitude, wavelength). I want to plot x, y, and amplitude on a 3D plot and have the color of the points correspond to the wavelength. I have seen many examples using functions to define color but my wavelength cannot be expressed by an analytic function. Is there a simple way to do this? Answer Here a another possible way to visualize 4D data: data = Flatten[Table[{x, y, x^2 + y^2, Sin[x - y]}, {x, -Pi, Pi,Pi/10}, {y,-Pi,Pi, Pi/10}], 1]; You can use the function Point along with VertexColors . Now the points are places using the first three elements and the color is determined by the fourth. In this case I used Hue, but you can use whatever you prefer. Graphics3D[ Point[data[[All, 1 ;; 3]], VertexColors -> Hue /@ data[[All, 4]]], Axes -> True, BoxRatios -> {1, 1, 1/GoldenRatio}]
Post a Comment
Read more

plotting - Filling between two spheres in SphericalPlot3D

Manipulate[ SphericalPlot3D[{1, 2 - n}, {θ, 0, Pi}, {ϕ, 0, 1.5 Pi}, Mesh -> None, PlotPoints -> 15, PlotRange -> {-2.2, 2.2}], {n, 0, 1}] I cant' seem to be able to make a filling between two spheres. I've already tried the obvious Filling -> {1 -> {2}} but Mathematica doesn't seem to like that option. Is there any easy way around this or ... Answer There is no built-in filling in SphericalPlot3D . One option is to use ParametricPlot3D to draw the surfaces between the two shells: Manipulate[ Show[SphericalPlot3D[{1, 2 - n}, {θ, 0, Pi}, {ϕ, 0, 1.5 Pi}, PlotPoints -> 15, PlotRange -> {-2.2, 2.2}], ParametricPlot3D[{ r {Sin[t] Cos[1.5 Pi], Sin[t] Sin[1.5 Pi], Cos[t]}, r {Sin[t] Cos[0 Pi], Sin[t] Sin[0 Pi], Cos[t]}}, {r, 1, 2 - n}, {t, 0, Pi}, PlotStyle -> Yellow, Mesh -> {2, 15}]], {n, 0, 1}]
Post a Comment
Read more

plotting - Mathematica: 3D plot based on combined 2D graphs

I have several sigmoidal fits to 3 different datasets, with mean fit predictions plus the 95% confidence limits (not symmetrical around the mean) and the actual data. I would now like to show these different 2D plots projected in 3D as in but then using proper perspective. In the link here they give some solutions to combine the plots using isometric perspective, but I would like to use proper 3 point perspective. Any thoughts? Also any way to show the mean points per time point for each series plus or minus the standard error on the mean would be cool too, either using points+vertical bars, or using spheres plus tubes. Below are some test data and the fit function I am using. Note that I am working on a logit(proportion) scale and that the final vertical scale is Log10(percentage). (* some test data *) data = Table[Null, {i, 4}]; data[[1]] = {{1, -5.8}, {2, -5.4}, {3, -0.8}, {4, -0.2}, {5, 4.6}, {1, -6.4}, {2, -5.6}, {3, -0.7}, {4, 0.04}, {5, 1.0}, {1, -6.8}, {2, -4.7}, {3, -1....
Post a Comment
Read more