0
$\begingroup$

I try to plot a phase diagram with ParametricPlot but everytime Mathematica does nothing.

I start by ;

sol = Solve[{α A k^(α - 1) - ρ - δ == 0,  A k^α - δ k - c == 0}, {c, k}] // FullSimplify


{csol, ksol} = {sol[[1, 1, 2]], sol[[1, 2, 2]]}


Clear[dec, dek, α, A, σ, ρ, δ];

dec = c'[t] == -(c[t]/σ) (A α k[t] ^(α - 1) - δ - ρ);

dek = k'[t] == -(A k[t]^α - c[t] - δ k[t]);

Clear[α, σ, ρ, δ];
α = 0.3; σ = 1; ρ = 0.04; δ = 0.02; A = 1;

Clear[nsba, cb, kb, cf, kf];
nsba = NDSolve[  {dec, dek, c[0] == csol (1 - 10^-5), k[0] == 
                                    ksol (1 - 10^-5)}, {c[t], k[t]}, {t, 0, 255}]

cb[t_] := Evaluate[nsba[[1, 1, 2]]]
kb[t_] := Evaluate[nsba[[1, 2, 2]]]

Clear[tende];
tende = 164;
cf[t_] := cb[tende - t]
kf[t_] := kb[tende - t]

After all, I try this last command ;

tra = ParametricPlot[{kf[t], cf[t]}, {t, 0, tende},
                     AxesLabel -> {k, c}, PlotStyle -> Thickness[0.01]];

Mathematica stays actionless. What could be the eventual problem ?

Improve this question
$\endgroup$
7
  • 1
    $\begingroup$ is the first Solve[ ] working for you? What is returning as solutions? $\endgroup$
    – Dr. belisarius
    Jul 21, 2015 at 22:32
  • 2
    $\begingroup$ The semicolon after your plot command suppresses output. $\endgroup$
    – Simon Rochester
    Jul 21, 2015 at 22:33
  • $\begingroup$ @belisarius yes the command Solve[ ] works fine. It gives me the right solution for k and c that I have also solved by hand. $\endgroup$
    – optimal control
    Jul 21, 2015 at 22:38
  • 1
    $\begingroup$ Please post the solutions here. Thanks! $\endgroup$
    – Dr. belisarius
    Jul 21, 2015 at 22:41
  • $\begingroup$ Optimal, unfortunately very little works in your code. Even trying to understand what you want to do by reading the code is pretty hard. Maybe you could describe the problem you are trying to solve in your question as well. As @belisarius hinted in his comment, the first Solve doesn't seem to work for me as written (it throws multiple errors). $\endgroup$
    – MarcoB
    Jul 21, 2015 at 23:20

1 Answer 1

Reset to default
2
$\begingroup$

There seem to be problems with your use of both Solve and NDSolve. However, ignoring the warnings Mathematica issues and plunging ahead, I get

sol = 
  Solve[{α A k^(α - 1) - ρ - δ == 0, A k^α - δ k - c == 0}, {c, k}] // FullSimplify;
{csol, ksol} = {sol[[1, 1, 2]], sol[[1, 2, 2]]};
dec = c'[t] == -(c[t]/σ) (A α k[t]^(α - 1) - δ - ρ);
dek = k'[t] == -(A k[t]^α - c[t] - δ k[t]);
α = 0.3; σ = 1; ρ = 0.04; δ = 0.02; A = 1;
nsba = 
  NDSolve[{dec, dek, c[0] == csol (1 - 10^-5), k[0] == ksol (1 - 10^-5)}, 
    {c[t], k[t]}, {t, 0, 255}];

At this point I arrive at the point in the code that is affecting your plot. I fix that by defining

cb = Head @ nsba[[1, 1, 2]];
kb = Head @ nsba[[1, 2, 2]];

Now

With[{tende = 164.},
  With[{cf = cb[tende - #] &, kf = kb[tende - #] &},
    ParametricPlot[{kf[t], cf[t]}, {t, 0, tende},
      AxesLabel -> {k, c},
      AspectRatio -> 1,
      PlotStyle -> Thickness[0.01]]]]

produces

plot

I have no idea whether or not this is a meaningful plot, but at least it shows you what your code is generating.

Improve this answer
$\endgroup$
1
  • $\begingroup$ Thanks so much, it is very helpful :) I think the problem is possible to be solved analytically. $\endgroup$
    – optimal control
    Jul 22, 2015 at 14:32

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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