       Re: Solving this in mathematica?

• To: mathgroup at smc.vnet.net
• Subject: [mg128066] Re: Solving this in mathematica?
• From: Bob Hanlon <hanlonr357 at gmail.com>
• Date: Thu, 13 Sep 2012 03:38:04 -0400 (EDT)
• Delivered-to: l-mathgroup@mail-archive0.wolfram.com
• Delivered-to: l-mathgroup@wolfram.com
• Delivered-to: mathgroup-newout@smc.vnet.net
• Delivered-to: mathgroup-newsend@smc.vnet.net
• References: <20120912065824.420E96824@smc.vnet.net>

```rbar = 0.006236 // Rationalize;

(* rt=r_bar You cannot use _ in a name.
I assume the RHS should be rbar *)

rt = rbar; (* rt is not used *)
k =
0.95 // Rationalize;
sigmar = 0.002 // Rationalize;
betazr = -0.00014 // Rationalize;
sigmaz = 0.4 // Rationalize;
pi = 0.99 // Rationalize; (* pi is not used *)
chi =
0.05 // Rationalize; (* chi is not used *)
Cbar = -3.7 // Rationalize;

(* Use RSolve to solve your recurrence equations *)

Clear[alpha1, alpha2];
sol1 = {alpha1[n], alpha2[n]} /. RSolve[{
alpha1[n] == alpha1[n - 1] + alpha2[n - 1],
alpha2[n] == k (alpha2[n - 1]),
alpha1 == 0, alpha2 == 1},
{alpha1[n], alpha2[n]}, n][];

alpha1[n_] = sol1[];
alpha2[n_] = sol1[];

Clear[sigma1sq, sigma12, sigma2sq];
sol2 = {sigma1sq[n], sigma12[n], sigma2sq[n]} /. RSolve[{
sigma1sq[n] == sigma2sq[n - 1] + 2 sigma12[n - 1] + sigmaz^2,
sigma12[n] == k (sigma12[n - 1]) + k (sigma2sq[n - 1]) + betazr,
sigma2sq[n] == (k^2) (sigma2sq[n - 1]) + sigmar^2,
sigma1sq == 0, sigma12 == 0, sigma2sq == 0},
{sigma1sq[n], sigma12[n], sigma2sq[n]}, n][] //
FullSimplify;

sigma1sq[n_] = sol2[];
sigma12[n_] = sol2[];
sigma2sq[n_] = sol2[];

(* You did not specify a starting value for psi[n], I have assumed that \
psi = 0 *)

Clear[phi1, phi2, psi];
sol3 = {phi1[n], phi2[n], psi[n]} /. RSolve[{
phi1[n] == phi1[n - 1] + phi2[n - 1] + (sigmaz^2)/2,
phi2[n] == k (phi2[n - 1]) + (1 - k) (rbar),
psi[n] == phi1[n] - sigma1sq[n]/2,
phi1 == 0, phi2 == 0, psi == 0},
{phi1[n], phi2[n], psi[n]}, n][] //
FullSimplify;

phi1[n_] = sol3[];
phi2[n_] = sol3[];
psi[n_] = sol3[];

B[h_, r_] = Exp[-alpha1[h]* r - psi[h]] //
PiecewiseExpand // Simplify;

NMinimize[Abs[Exp[Cbar - beta] Sum[(Pi^x) B[x, r],
{x, 1, 1000}] - 1], {r, beta},
WorkingPrecision -> 35]

{1.0000000000000000000000000000000000, {r ->
311.81988579682125561865976837161440,
beta -> 158.37558965184462887837320791605232}}

This apears to indicate that there is no solution to your equation
since the minimum is not zero. Although this may be just a local
minimum and the function may need further exploration.

Bob Hanlon

On Wed, Sep 12, 2012 at 2:58 AM, Leon <leon.he88 at gmail.com> wrote:
> I have following code in mathematica:
>
> ------------------------------------------------------------
> rbar = 0.006236
> rt = r_bar
> k = 0.95
> sigmar = 0.002
> betazr = -0.00014
> sigmaz = 0.4
> pi = 0.99
> chi = 0.05
> Cbar = -3.7
>
> alpha1[n_] := alpha1[n] = alpha1[n - 1] + alpha2[n - 1]
> alpha2[n_] := alpha2[n] = k (alpha2[n - 1])
> sigma1sq[n_] :=
>  sigma1sq[n] = sigma2sq[n - 1] + 2 sigma12[n - 1] + sigmaz^2
> sigma12[n_] :=
>  sigma12[n] = k (sigma12[n - 1]) + k (sigma2sq[n - 1]) + betazr
> sigma2sq[n_] := sigma2sq[n] = (k^2) (sigma2sq[n - 1]) + sigmar^2
> phi1[n_] := phi1[n] = phi1[n - 1] + phi2[n - 1] + (0.5) (sigmaz^2)
> phi2[n_] := phi2[n] = k (phi2[n - 1]) + (1 - k) (rbar)
> psi[n_] := psi[n] = phi1[n] - (0.5) (sigma1sq[n])
>
> alpha1 = 0
> alpha2 = 1
> sigma1sq = 0
> sigma12 = 0
> sigma2sq = 0
> phi1 = 0
> phi2 = 0
>
> B[h_, r_] := Exp[(-alpha1[h]) (r) - psi[h]]
> Exp[Cbar - beta] Sum[(Pi^x) B[x, r], {x, 1, 1000}]
>
> ----------------------------------------------------------------------
>
> and I am wondering if it is possible to solve the last line such that I have "r" as a function of "beta", satisfying
>
> Exp[Cbar - beta] Sum[(Pi^x) B[x, r], {x, 1, 1000}] = 1
>
> Because ultimately, I would need to integrate a function J[r] over "beta", so if I don't have "r" as a function of "beta", I don't know how to do the integration of J[r].
>
> Thank you!! L.
>

```

• Prev by Date: Eigenvalue and eigenvectors of a 10x10 matrix
• Next by Date: Re: Are some equations unsolvable?
• Previous by thread: Solving this in mathematica?
• Next by thread: Re: Solving this in mathematica?