1
$\begingroup$

I've tried to draw a vector plot with piecewise functions, but it shows empty result.
When I don't use piecewise functions, Vectorplot works well..
I can't figure out what the problem is.
I would really appreciate if someone helps me.
The whole code is below.

A = 5;
b = 0.6;
theta = 0.7;
lambda = 1;
rho = 0.9;
eta = 0.3;
B = 2;
alpha = 1;
i = 0.1;
uppergamma = 1.4;
lowergamma = 0.4;
Delta = 1;
delta = 1;

u[q_] = A q^b;
up[q_] = A b q^(b-1);

F[x_] = 
Piecewise[{
{1, x > uppergamma},
{0, x < lowergamma},
{(x-lowergamma)/(uppergamma-lowergamma), lowergamma <= x <= uppergamma}
}];

integrateF[x_] = 
Piecewise[{
{0, x < lowergamma},
{x-uppergamma + 0.5, x > uppergamma},
{0.5*(x-lowergamma)^2, lowergamma <= x <= uppergamma}
}];

yinner[n_, i_] = y /. FindRoot[theta*(up[y]-1)/((1-theta)*up[y]+theta) == i/(alpha*n), {y, 0.001}];
y[n_,i_] = If[Im[yinner[n,i]]==0, yinner[n,i], 0];

niso[gamma_] = Max[0, lambda*F[gamma]/(delta+lambda*F[gamma])];
gammaiso[n_] = gamma /. 
FindRoot[(rho+delta)*gamma + lambda*integrateF[gamma] 
== alpha*(1-theta)*(u[y[n,i]]-y[n,i]), {gamma, 0.001}]

Show[{


ParametricPlot[{n, gammaiso[n]}, {n, 0, 1},  PlotRange -> {{0,1}, {0,uppergamma}},
PlotStyle -> {Blue}, AxesLabel -> {n, gamma},
BaseStyle -> {FontSize -> 12}, PlotLegends -> {Style["gamma_R dot = 0", 12]}],


ParametricPlot[{niso[gamma], gamma}, {gamma, 0, uppergamma}, PlotRange -> {{0,1}, {0,uppergamma}},
PlotStyle -> {Red}, AxesLabel -> {n, gamma},
BaseStyle -> {FontSize -> 12}, PlotLegends -> {Style["ndot=0", 12]}], 


VectorPlot[{(1-n)*lambda*F[gamma]-delta*n,
(rho+delta)*gamma + lambda*integrateF[gamma]  - alpha *(1-theta)*(u[y[n,i]]-y[n,i])},
{n, 0, 1}, {gamma, 0, uppergamma}]

}]

enter image description here

$\endgroup$
2
  • 3
    $\begingroup$ I just replaced all your = with := and it worked. Why are you defining functions with = and not :=? always use := for functions unless you have specific reason not to. V 13.1 screen shot !Mathematica graphics $\endgroup$
    – Nasser
    Jul 24 at 6:00
  • 3
    $\begingroup$ Just to add to Nasser’s useful comment: it’s ok—preferred even—to use = on your constants (as in b=6). As Nasser says, use := in function defs (e.g., f[x_]:=x). There are exceptions that have to do with efficiency, but for most users := is the way to go for function definitions. Trying to understand why would be a rewarding task. $\endgroup$ Jul 24 at 6:13

1 Answer 1

0
$\begingroup$
Clear["Global`*"]

Use exact values for constants.

A = 5;
b = 3/5;
theta = 7/10;
lambda = 1;
rho = 9/10;
eta = 3/10;
B = 2;
alpha = 1;
i = 1/10;
uppergamma = 7/5;
lowergamma = 2/5;
Delta = 1;
delta = 1;

u[q_] = A q^b;
up[q_] = A b q^(b - 1);

F[x_] = Piecewise[{{1, x > uppergamma}, {0, 
     x < lowergamma}, {(x - lowergamma)/(uppergamma - lowergamma), 
     lowergamma <= x <= uppergamma}}];

integrateF[x_] = 
  Piecewise[{{0, x < lowergamma}, {x - uppergamma + 1/2, 
      x > uppergamma}, {1/2*(x - lowergamma)^2, 
      lowergamma <= x <= uppergamma}}] // Simplify;

yinner can be solved exactly

yinner[n_, iv_] = 
 SolveValues[theta*(up[y] - 1)/((1 - theta)*up[y] + theta) == iv/(alpha*n), y,
     Reals][[1]] // Simplify

enter image description here

The condition on yinner indicates when it is real

y[n_, iv_] = Piecewise[{List @@ yinner[n, iv]}]

enter image description here

niso[gamma_] = 
  Max[0, lambda*F[gamma]/(delta + lambda*F[gamma])] // PiecewiseExpand // 
   FullSimplify;

Functions such as FindRoot that use numeric techniques should restrict their arguments to numeric values and use SetDelayed.

gammaiso[n_?NumericQ] := 
 gamma /. FindRoot[(rho + delta)*gamma + lambda*integrateF[gamma] == 
    alpha*(1 - theta)*(u[y[n, i]] - y[n, i]), {gamma, 1/1000}, 
   WorkingPrecision -> 15]

Plotting,

Show[{ParametricPlot[{n, gammaiso[n]}, {n, 0, 1}, 
   PlotRange -> {{0, 1}, {0, uppergamma}}, PlotStyle -> {Blue}, 
   AxesLabel -> {n, gamma}, BaseStyle -> {FontSize -> 12}, 
   PlotLegends -> {Style["gamma_R dot = 0", 12]}], 
  ParametricPlot[{niso[gamma], gamma}, {gamma, 0, uppergamma}, 
   PlotRange -> {{0, 1}, {0, uppergamma}}, PlotStyle -> {Red}, 
   AxesLabel -> {n, gamma}, BaseStyle -> {FontSize -> 12}, 
   PlotLegends -> {Style["ndot=0", 12]}], 
  VectorPlot[{(1 - n)*lambda*F[gamma] - delta*n, (rho + delta)*gamma + 
        lambda*integrateF[gamma] - alpha*(1 - theta)*(u[y[n, i]] - y[n, i])} //
       PiecewiseExpand // FullSimplify // Evaluate, {n, 0, 1}, {gamma, 0, 
    uppergamma}]}]

enter image description here

$\endgroup$

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.