Advanced Search

Maple Questions and Posts

These are Posts and Questions associated with the product, Maple
View unanswered posts
Maple Questions and Posts Feed

Kernel connection lost...

Asked by: jensw87 5
April 08 2015
1  4

Hello,

I'm trying to calculate and export a huge "Strategy-Matrix". The calculation works in a smaller version with less possible strategies. For this version I use a specialised simulation-PC with enough memory (that was my problem using my own PC). But now I've got this error massage ("Kernel connection has been lost. You should save this worksheet and restart Maple. Executing commands in Maple requires a connection to the Maple kernel. Firewalls have been known to cause problems with kernel connections in Maple. Please ensure that any firewall software is configured to allow Maple to create connections to the kernel. Consult the FAQ for more information.")

The firewall should not be the problem. Below the code. num_Strat is 63001. 

Thank you in advance for your help.

 

m_Gi:= Matrix(num_Strat, num_Strat):

for t from 1 by 1 to num_Strat
do
Digits:= 5:
v_h2:= evalf((-(24/17)* v_pEi[t]^2)*v_h1 + (v_pEi[t]*(-(128/51)*v_pKi[t]+980/51))*v_h1+ v_pEi[t]*((40/51)*v_pKj+(15/17)*v_pEj)-((110/51)*v_pKi[t]^2)*v_h1 + v_pKi[t]*((55/51)*v_pKj+1220/51*v_h1+(10/17)*v_pEj)):
m_Gi(..,t):= v_h2:
end do:
ExportMatrix(TestMatrixNK, m_Gi):

plot this series...

Asked by: Parham2016 70
plot series
April 08 2015
0  1

plot.mw

 

Hi all. Please help me to plot this. I don't know why it doesn't work correctly.

Regards

pdsolve/Solve for derivatives algebraically (witho...

Asked by: jschulzb 150
pde pdsolve
April 07 2015
0  1

Hi, 

     I have a question regarding pdsolve, or Solve from the PDEtools package. I have a set of equations relating partial derivatives, and I'd like to isolate certain terms without explicitly known the functions. I can do this for a single equation, but not multiple ones. I'm curious if Maple can currently handle a system of eqns like these easily, since I will be increasing the number of eqns in the future. Here's the code 

 

 

> 

restart;
> 

PDEtools:-declare(H=H(x,y,t)):

H(x, y, t)*`will now be displayed as`*H

 (1)
> 

eq1:= H[tt](x,y,t) = H[xx](x,y,t) + H[yy](x,y,t);

H[tt](x, y, t) = H[xx](x, y, t)+H[yy](x, y, t)

 (2)
> 

eq2 := diff(H[tt](x,y,t), t) = diff(H[tx](x,y,t), x) + diff(H[ty](x,y,t), y);

diff(H[tt](x, y, t), t) = diff(H[tx](x, y, t), x)+diff(H[ty](x, y, t), y)

 (3)
> 

eq3 := diff(H[tx](x,y,t), t) = diff(H[xx](x,y,t), x) + diff(H[xy](x,y,t), y);

diff(H[tx](x, y, t), t) = diff(H[xx](x, y, t), x)+diff(H[xy](x, y, t), y)

 (4)
> 

eq4 :=diff(H[ty](x,y,t), t) = diff(H[xy](x,y,t), x) + diff(H[yy](x,y,t), y);

diff(H[ty](x, y, t), t) = diff(H[xy](x, y, t), x)+diff(H[yy](x, y, t), y)

 (5)
> 

PDEtools:-Solve(eq3, H[xy]);

H[xy](x, y, t) = Int(diff(H[tx](x, y, t), t)-(diff(H[xx](x, y, t), x)), y)+_F1(x, t)

 (6)
> 

PDEtools:-Solve({eq1, eq2, eq3, eq4}, H[xy]);

Error, (in pdsolve/sys) the input system cannot contain equations in the arbitrary parameters alone; found equation depending only on {H[tt](x,y,t), H[xx](x,y,t), H[yy](x,y,t)}: H[tt](x,y,t)-H[xx](x,y,t)-H[yy](x,y,t)

  
> 

 
> 

 

 

Download PDESolveHelp.mw

Combinations of numerical solutions from dsolve...

Asked by: bogo 55
dsolve stiff
April 07 2015
0  3

Hi! I am trying to plot and store in memory some specific combinations of the solutions of the systems of ODEs that I get numerically from dsolve for a particular range of the independent variable. 

A particular case for my problem is the following system of stiff ODEs for two unknown functions f[0,0](x) and f[1,0](x) beween xini (where the Initial conditions are defined) and xfin, an arbitrary value of x. Note that rosebrock method does not work, and I can only solve it with lsode[adamsfull] or lsode[backfull]. I am attaching a maple file that shows what I have done.


``

############## System of ODEs that needs to be solved ####################################

xini := .1

.1

 (1)

xfin := 2

2

 (2)

SystemToSolve := diff(f[0, 0](x), x)+(2./x^5+.5000000000/x)*f[0, 0](x) = -15.58845727*sin(.5773502693*x)/x^2+(46.76537182*(3.*sin(.5773502693*x)-1.732050808*x*cos(.5773502693*x)))/x^4, diff(f[1, 0](x), x)+(6./x^5+1.500000000/x)*f[1, 0](x)-1.*f[0, 0](x)/x = (-15.58845727*sin(.5773502693*x)/x^2+(46.76537182*(3.*sin(.5773502693*x)-1.732050808*x*cos(.5773502693*x)))/x^4)*(1.-1.*(1/x^4)^(1/4)*exp(1/x^4)*GAMMA(.7500000000, 1/x^4))-(1.*(-10.39230485*sin(.5773502693*x)/x^2+(31.17691454*(3.*sin(.5773502693*x)-1.732050808*x*cos(.5773502693*x)))/x^4+(4.*((.8660254040*(3.*sin(.5773502693*x)-1.732050808*x*cos(.5773502693*x)))/x+(.8660254040*((3.*(1.-6./x^2))*sin(.5773502693*x)+10.39230485*cos(.5773502693*x)/x))/x))/((1/x^4)^(1/4)*exp(1/x^4)*GAMMA(.7500000000, 1/x^4)*x^5)))*(1/x^4)^(1/4)*exp(1/x^4)*GAMMA(.7500000000, 1/x^4), f[0, 0](xini) = 1.503498546, f[1, 0](xini) = -.5011661819:

 

###################################################################################

 

``

ListProcs := dsolve({SystemToSolve}, numeric, method = lsode[backfull], output = listprocedure):

f00 := eval(f[0, 0](x), ListProcs);

proc (x) local _res, _dat, _solnproc, _xout, _ndsol, _pars, _i; option `Copyright (c) 2000 by Waterloo Maple Inc. All rights reserved.`; if 1 < nargs then error "invalid input: too many arguments" end if; _EnvDSNumericSaveDigits := Digits; Digits := 15; if _EnvInFsolve = true then _xout := evalf[_EnvDSNumericSaveDigits](x) else _xout := evalf(x) end if; _dat := Array(1..4, {(1) = proc (_xin) local _xout, _n, _y0, _ctl, _octl, _reinit, _errcd, _fcn, _i, _yini, _pars, _ini, _par; option `Copyright (c) 2002 by the University of Waterloo. All rights reserved.`; table( [( "complex" ) = false ] ) _xout := _xin; _ctl := array( 1 .. 34, [( 1 ) = (2), ( 2 ) = (.1), ( 3 ) = (.1), ( 4 ) = (1), ( 5 ) = (1), ( 6 ) = (22), ( 7 ) = (0), ( 9 ) = (-.5011661819), ( 8 ) = (1.503498546), ( 11 ) = (0.1e-6), ( 10 ) = (0.1e-6), ( 13 ) = (0), ( 12 ) = (0), ( 15 ) = (0), ( 14 ) = (0), ( 18 ) = (0), ( 19 ) = (0), ( 16 ) = (0), ( 17 ) = (0), ( 22 ) = (0), ( 23 ) = (0), ( 20 ) = (0), ( 21 ) = (0), ( 27 ) = (0), ( 26 ) = (0), ( 25 ) = (0), ( 24 ) = (0), ( 31 ) = (0), ( 30 ) = (0), ( 29 ) = (0), ( 28 ) = (0), ( 32 ) = (0), ( 33 ) = (-1), ( 34 ) = (0) ] ); _octl := array( 1 .. 34, [( 1 ) = (2), ( 2 ) = (.1), ( 3 ) = (.1), ( 4 ) = (1), ( 5 ) = (1), ( 6 ) = (22), ( 7 ) = (0), ( 9 ) = (-.5011661819), ( 8 ) = (1.503498546), ( 11 ) = (0.1e-6), ( 10 ) = (0.1e-6), ( 13 ) = (0), ( 12 ) = (0), ( 15 ) = (0), ( 14 ) = (0), ( 18 ) = (0), ( 19 ) = (0), ( 16 ) = (0), ( 17 ) = (0), ( 22 ) = (0), ( 23 ) = (0), ( 20 ) = (0), ( 21 ) = (0), ( 27 ) = (0), ( 26 ) = (0), ( 25 ) = (0), ( 24 ) = (0), ( 31 ) = (0), ( 30 ) = (0), ( 29 ) = (0), ( 28 ) = (0), ( 32 ) = (0), ( 33 ) = (-1), ( 34 ) = (0) ] ); _n := trunc(_ctl[1]); _yini := Array(0..2, {(1) = .1, (2) = 1.503498546}); _y0 := Array(0..2, {(1) = .1, (2) = 1.503498546}); _fcn := proc (N, X, Y, YP) option `[Y[1] = f[0,0](x), Y[2] = f[1,0](x)]`; YP[1] := -15.58845727*sin(.5773502693*X)/X^2+46.76537182*(3.*sin(.5773502693*X)-1.732050808*X*cos(.5773502693*X))/X^4-(2./X^5+.5000000000/X)*Y[1]; if 1/X^4 < 0 then YP[1] := undefined; return 0 end if; if 1/X^4 < 0 then YP[1] := undefined; return 0 end if; YP[2] := (-15.58845727*sin(.5773502693*X)/X^2+46.76537182*(3.*sin(.5773502693*X)-1.732050808*X*cos(.5773502693*X))/X^4)*(1.-1.*evalf((1/X^4)^(1/4))*exp(1/X^4)*GAMMA(.7500000000, 1/X^4))-1.*(-10.39230485*sin(.5773502693*X)/X^2+31.17691454*(3.*sin(.5773502693*X)-1.732050808*X*cos(.5773502693*X))/X^4+4.*(.8660254040*(3.*sin(.5773502693*X)-1.732050808*X*cos(.5773502693*X))/X+.8660254040*(3.*(1.-6./X^2)*sin(.5773502693*X)+10.39230485*cos(.5773502693*X)/X)/X)*evalf(1/(1/X^4)^(1/4))/(exp(1/X^4)*GAMMA(.7500000000, 1/X^4)*X^5))*evalf((1/X^4)^(1/4))*exp(1/X^4)*GAMMA(.7500000000, 1/X^4)-(6./X^5+1.500000000/X)*Y[2]+1.*Y[1]/X; 0 end proc; _pars := []; if not type(_xout, 'numeric') then if member(_xout, ["start", "left", "right"]) then return _y0[0] elif _xout = "method" then return "lsode" elif _xout = "numfun" then return trunc(_ctl[24+trunc(_ctl[1])]) elif _xout = "initial" then return [seq(_yini[_i], _i = 0 .. _n)] elif _xout = "parameters" then return [seq(_yini[_n+_i], _i = 1 .. nops(_pars))] elif _xout = "initial_and_parameters" then return [seq(_yini[_i], _i = 0 .. _n)], [seq(_yini[_n+_i], _i = 1 .. nops(_pars))] elif _xout = "last" then if _ctl[2]-_y0[0] = 0. then error "no information is available on last computed point" else _xout := _ctl[2] end if elif _xout = "enginedata" then return eval(_octl, 1) elif _xout = "function" then return eval(_fcn, 1) elif type(_xin, `=`) and type(rhs(_xin), 'list') and member(lhs(_xin), {"initial", "parameters", "initial_and_parameters"}) then _ini, _par := [], []; if lhs(_xin) = "initial" then _ini := rhs(_xin) elif lhs(_xin) = "parameters" then _par := rhs(_xin) elif select(type, rhs(_xin), `=`) <> [] then _par, _ini := selectremove(type, rhs(_xin), `=`) elif nops(rhs(_xin)) < nops(_pars)+1 then error "insufficient data for specification of initial and parameters" else _par := rhs(_xin)[-nops(_pars) .. -1]; _ini := rhs(_xin)[1 .. -nops(_pars)-1] end if; _xout := lhs(_xout); if _par <> [] then `dsolve/numeric/process_parameters`(_n, _pars, _par, _yini) end if; if _ini <> [] then `dsolve/numeric/process_initial`(_n, _ini, _yini, _pars) end if; if _pars <> [] then _par := {seq(rhs(_pars[_i]) = _yini[_n+_i], _i = 1 .. nops(_pars))}; for _i from 0 to _n do _y0[_i] := subs(_par, _yini[_i]) end do; for _i from _n+1 to _n+nops(_pars) do _y0[_i] := _yini[_i] end do else for _i from 0 to _n do _y0[_i] := _yini[_i] end do end if; _octl[2] := _y0[0]; _octl[3] := _y0[0]; for _i to _n do _octl[_i+7] := _y0[_i] end do; for _i to nops(_pars) do _octl[2*_n+30+_i] := _y0[_n+_i] end do; for _i to 34 do _ctl[_i] := _octl[_i] end do; if _Env_smart_dsolve_numeric = true and type(_y0[0], 'numeric') then procname("right") := _y0[0]; procname("left") := _y0[0] end if; if _xout = "initial" then return [seq(_yini[_i], _i = 0 .. _n)] elif _xout = "parameters" then return [seq(_yini[_n+_i], _i = 1 .. nops(_pars))] else return [seq(_yini[_i], _i = 0 .. _n)], [seq(_yini[_n+_i], _i = 1 .. nops(_pars))] end if else return "procname" end if end if; if _xout-_y0[0] = 0. then return [seq(_y0[_i], _i = 0 .. _n)] end if; _reinit := false; if _xin <> "last" then if 0 < 0 and `dsolve/numeric/checkglobals`(0, table( [ ] ), _pars, _n, _yini) then _reinit := true; if _pars <> [] then _par := {seq(rhs(_pars[_i]) = _yini[_n+_i], _i = 1 .. nops(_pars))}; for _i from 0 to _n do _y0[_i] := subs(_par, _yini[_i]) end do; for _i from _n+1 to _n+nops(_pars) do _y0[_i] := _yini[_i] end do else for _i from 0 to _n do _y0[_i] := _yini[_i] end do end if; for _i to _n do _octl[_i+7] := _y0[_i] end do; for _i to nops(_pars) do _octl[2*_n+30+_i] := _y0[_n+_i] end do end if; if _pars <> [] and select(type, {seq(_yini[_n+_i], _i = 1 .. nops(_pars))}, 'undefined') <> {} then error "parameters must be initialized before solution can be computed" end if end if; if not _reinit and _xout-_ctl[2] = 0 then [_ctl[2], seq(_ctl[_i], _i = 8 .. 7+_n)] else if sign(_xout-_ctl[2]) <> sign(_ctl[2]-_y0[0]) or abs(_xout-_y0[0]) < abs(_xout-_ctl[2]) or _reinit then for _i to 34 do _ctl[_i] := _octl[_i] end do end if; _ctl[3] := _xout; if Digits <= evalhf(Digits) then try _errcd := evalhf(`dsolve/numeric/lsode`(_fcn, var(_ctl))) catch: userinfo(2, `dsolve/debug`, print(`Exception in lsode:`, [lastexception])); if searchtext('evalhf', lastexception[2]) <> 0 or searchtext('real', lastexception[2]) <> 0 or searchtext('hardware', lastexception[2]) <> 0 then _errcd := `dsolve/numeric/lsode`(_fcn, _ctl) else error end if end try else _errcd := `dsolve/numeric/lsode`(_fcn, _ctl) end if; if _errcd < 0 then userinfo(2, {dsolve, `dsolve/lsode`}, `Last values returned:`); userinfo(2, {dsolve, `dsolve/lsode`}, ` t =`, _ctl[2]); _i := 8; userinfo(2, {dsolve, `dsolve/lsode`}, ` y =`, _ctl[_i]); for _i from _i+1 to 7+_n do userinfo(2, {dsolve, `dsolve/lsode`}, ` `, _ctl[_i]) end do; if _errcd+1. = 0. then if _ctl[14+trunc(_ctl[1])] <> 0 then error "an excessive amount of work was done, maxstep may be too small" else error "an excessive amount of work (greater than mxstep) was done" end if elif _errcd+2. = 0. then error "too much accuracy was requested for the machine being used" elif _errcd+3. = 0. then error "illegal input was detected" elif _errcd+4. = 0. then error "repeated error test failures on the attempted step" elif _errcd+5. = 0. then error "repeated convergence test failures on the attempted step" elif _errcd+6. = 0. then error "pure relative error control requested for a variable that has vanished" elif _errcd+7. = 0. then error "cannot evaluate the solution past %1, maxfun limit exceeded (see <a href='http://www.maplesoft.com/support/help/search.aspx?term=dsolve,maxfun' target='_new'>?dsolve,maxfun</a> for details)", evalf[8](_ctl[2]) else error "unknown error code returned from lsode %1", trunc(_errcd) end if end if; if _Env_smart_dsolve_numeric = true then if _y0[0] < _xout and procname("right") < _xout then procname("right") := _xout elif _xout < _y0[0] and _xout < procname("left") then procname("left") := _xout end if end if; [_xout, seq(_ctl[_i], _i = 8 .. 7+_n)] end if end proc, (2) = Array(1..3, {(1) = 18446744078356217278, (2) = 18446744078356217454, (3) = 18446744078356217630}), (3) = [x, f[0, 0](x), f[1, 0](x)], (4) = []}); _solnproc := _dat[1]; _pars := map(rhs, _dat[4]); if not type(_xout, 'numeric') then if member(x, ["start", 'start', "method", 'method', "left", 'left', "right", 'right', "leftdata", "rightdata", "enginedata", "eventstop", 'eventstop', "eventclear", 'eventclear', "eventstatus", 'eventstatus', "eventcount", 'eventcount', "laxtol", 'laxtol', "numfun", 'numfun', NULL]) then _res := _solnproc(convert(x, 'string')); if 1 < nops([_res]) then return _res elif type(_res, 'array') then return eval(_res, 1) elif _res <> "procname" then return _res end if elif member(x, ["last", 'last', "initial", 'initial', NULL]) then _res := _solnproc(convert(x, 'string')); if type(_res, 'list') then return _res[2] else return NULL end if elif member(x, ["parameters", 'parameters', "initial_and_parameters", 'initial_and_parameters', NULL]) then _xout := convert(x, 'string'); _res := _solnproc(_xout); if _xout = "parameters" then return [seq(_pars[_i] = _res[_i], _i = 1 .. nops(_pars))] else return [_res[2], seq(_pars[_i] = [_res][2][_i], _i = 1 .. nops(_pars))] end if elif type(_xout, `=`) and member(lhs(_xout), ["initial", 'initial', "parameters", 'parameters', "initial_and_parameters", 'initial_and_parameters', NULL]) then _xout := convert(lhs(x), 'string') = rhs(x); if lhs(_xout) = "initial" then if type(rhs(_xout), 'list') then _res := _solnproc(_xout) else _res := _solnproc("initial" = ["single", 2, rhs(_xout)]) end if elif not type(rhs(_xout), 'list') then error "initial and/or parameter values must be specified in a list" elif lhs(_xout) = "initial_and_parameters" and nops(rhs(_xout)) = nops(_pars)+1 then _res := _solnproc(lhs(_xout) = ["single", 2, op(rhs(_xout))]) else _res := _solnproc(_xout) end if; if lhs(_xout) = "initial" then return _res[2] elif lhs(_xout) = "parameters" then return [seq(_pars[_i] = _res[_i], _i = 1 .. nops(_pars))] else return [_res[2], seq(_pars[_i] = [_res][2][_i], _i = 1 .. nops(_pars))] end if elif type(_xout, `=`) and member(lhs(_xout), ["eventdisable", 'eventdisable', "eventenable", 'eventenable', "eventfired", 'eventfired', "direction", 'direction', NULL]) then return _solnproc(convert(lhs(x), 'string') = rhs(x)) elif _xout = "solnprocedure" then return eval(_solnproc) elif _xout = "sysvars" then return _dat[3] end if; if procname <> unknown then return ('procname')(x) else _ndsol := `tools/gensym`("f[0,0](x)"); eval(FromInert(_Inert_FUNCTION(_Inert_NAME("assign"), _Inert_EXPSEQ(ToInert(_ndsol), _Inert_VERBATIM(pointto(_dat[2][2])))))); return FromInert(_Inert_FUNCTION(ToInert(_ndsol), _Inert_EXPSEQ(ToInert(x)))) end if end if; try _res := _solnproc(_xout); _res[2] catch: error end try end proc

 (3)

f10 := eval(f[1, 0](x), ListProcs);

proc (x) local _res, _dat, _solnproc, _xout, _ndsol, _pars, _i; option `Copyright (c) 2000 by Waterloo Maple Inc. All rights reserved.`; if 1 < nargs then error "invalid input: too many arguments" end if; _EnvDSNumericSaveDigits := Digits; Digits := 15; if _EnvInFsolve = true then _xout := evalf[_EnvDSNumericSaveDigits](x) else _xout := evalf(x) end if; _dat := Array(1..4, {(1) = proc (_xin) local _xout, _n, _y0, _ctl, _octl, _reinit, _errcd, _fcn, _i, _yini, _pars, _ini, _par; option `Copyright (c) 2002 by the University of Waterloo. All rights reserved.`; table( [( "complex" ) = false ] ) _xout := _xin; _ctl := array( 1 .. 34, [( 1 ) = (2), ( 2 ) = (.1), ( 3 ) = (.1), ( 4 ) = (1), ( 5 ) = (1), ( 6 ) = (22), ( 7 ) = (0), ( 9 ) = (-.5011661819), ( 8 ) = (1.503498546), ( 11 ) = (0.1e-6), ( 10 ) = (0.1e-6), ( 13 ) = (0), ( 12 ) = (0), ( 15 ) = (0), ( 14 ) = (0), ( 18 ) = (0), ( 19 ) = (0), ( 16 ) = (0), ( 17 ) = (0), ( 22 ) = (0), ( 23 ) = (0), ( 20 ) = (0), ( 21 ) = (0), ( 27 ) = (0), ( 26 ) = (0), ( 25 ) = (0), ( 24 ) = (0), ( 31 ) = (0), ( 30 ) = (0), ( 29 ) = (0), ( 28 ) = (0), ( 32 ) = (0), ( 33 ) = (-1), ( 34 ) = (0) ] ); _octl := array( 1 .. 34, [( 1 ) = (2), ( 2 ) = (.1), ( 3 ) = (.1), ( 4 ) = (1), ( 5 ) = (1), ( 6 ) = (22), ( 7 ) = (0), ( 9 ) = (-.5011661819), ( 8 ) = (1.503498546), ( 11 ) = (0.1e-6), ( 10 ) = (0.1e-6), ( 13 ) = (0), ( 12 ) = (0), ( 15 ) = (0), ( 14 ) = (0), ( 18 ) = (0), ( 19 ) = (0), ( 16 ) = (0), ( 17 ) = (0), ( 22 ) = (0), ( 23 ) = (0), ( 20 ) = (0), ( 21 ) = (0), ( 27 ) = (0), ( 26 ) = (0), ( 25 ) = (0), ( 24 ) = (0), ( 31 ) = (0), ( 30 ) = (0), ( 29 ) = (0), ( 28 ) = (0), ( 32 ) = (0), ( 33 ) = (-1), ( 34 ) = (0) ] ); _n := trunc(_ctl[1]); _yini := Array(0..2, {(1) = .1, (2) = 1.503498546}); _y0 := Array(0..2, {(1) = .1, (2) = 1.503498546}); _fcn := proc (N, X, Y, YP) option `[Y[1] = f[0,0](x), Y[2] = f[1,0](x)]`; YP[1] := -15.58845727*sin(.5773502693*X)/X^2+46.76537182*(3.*sin(.5773502693*X)-1.732050808*X*cos(.5773502693*X))/X^4-(2./X^5+.5000000000/X)*Y[1]; if 1/X^4 < 0 then YP[1] := undefined; return 0 end if; if 1/X^4 < 0 then YP[1] := undefined; return 0 end if; YP[2] := (-15.58845727*sin(.5773502693*X)/X^2+46.76537182*(3.*sin(.5773502693*X)-1.732050808*X*cos(.5773502693*X))/X^4)*(1.-1.*evalf((1/X^4)^(1/4))*exp(1/X^4)*GAMMA(.7500000000, 1/X^4))-1.*(-10.39230485*sin(.5773502693*X)/X^2+31.17691454*(3.*sin(.5773502693*X)-1.732050808*X*cos(.5773502693*X))/X^4+4.*(.8660254040*(3.*sin(.5773502693*X)-1.732050808*X*cos(.5773502693*X))/X+.8660254040*(3.*(1.-6./X^2)*sin(.5773502693*X)+10.39230485*cos(.5773502693*X)/X)/X)*evalf(1/(1/X^4)^(1/4))/(exp(1/X^4)*GAMMA(.7500000000, 1/X^4)*X^5))*evalf((1/X^4)^(1/4))*exp(1/X^4)*GAMMA(.7500000000, 1/X^4)-(6./X^5+1.500000000/X)*Y[2]+1.*Y[1]/X; 0 end proc; _pars := []; if not type(_xout, 'numeric') then if member(_xout, ["start", "left", "right"]) then return _y0[0] elif _xout = "method" then return "lsode" elif _xout = "numfun" then return trunc(_ctl[24+trunc(_ctl[1])]) elif _xout = "initial" then return [seq(_yini[_i], _i = 0 .. _n)] elif _xout = "parameters" then return [seq(_yini[_n+_i], _i = 1 .. nops(_pars))] elif _xout = "initial_and_parameters" then return [seq(_yini[_i], _i = 0 .. _n)], [seq(_yini[_n+_i], _i = 1 .. nops(_pars))] elif _xout = "last" then if _ctl[2]-_y0[0] = 0. then error "no information is available on last computed point" else _xout := _ctl[2] end if elif _xout = "enginedata" then return eval(_octl, 1) elif _xout = "function" then return eval(_fcn, 1) elif type(_xin, `=`) and type(rhs(_xin), 'list') and member(lhs(_xin), {"initial", "parameters", "initial_and_parameters"}) then _ini, _par := [], []; if lhs(_xin) = "initial" then _ini := rhs(_xin) elif lhs(_xin) = "parameters" then _par := rhs(_xin) elif select(type, rhs(_xin), `=`) <> [] then _par, _ini := selectremove(type, rhs(_xin), `=`) elif nops(rhs(_xin)) < nops(_pars)+1 then error "insufficient data for specification of initial and parameters" else _par := rhs(_xin)[-nops(_pars) .. -1]; _ini := rhs(_xin)[1 .. -nops(_pars)-1] end if; _xout := lhs(_xout); if _par <> [] then `dsolve/numeric/process_parameters`(_n, _pars, _par, _yini) end if; if _ini <> [] then `dsolve/numeric/process_initial`(_n, _ini, _yini, _pars) end if; if _pars <> [] then _par := {seq(rhs(_pars[_i]) = _yini[_n+_i], _i = 1 .. nops(_pars))}; for _i from 0 to _n do _y0[_i] := subs(_par, _yini[_i]) end do; for _i from _n+1 to _n+nops(_pars) do _y0[_i] := _yini[_i] end do else for _i from 0 to _n do _y0[_i] := _yini[_i] end do end if; _octl[2] := _y0[0]; _octl[3] := _y0[0]; for _i to _n do _octl[_i+7] := _y0[_i] end do; for _i to nops(_pars) do _octl[2*_n+30+_i] := _y0[_n+_i] end do; for _i to 34 do _ctl[_i] := _octl[_i] end do; if _Env_smart_dsolve_numeric = true and type(_y0[0], 'numeric') then procname("right") := _y0[0]; procname("left") := _y0[0] end if; if _xout = "initial" then return [seq(_yini[_i], _i = 0 .. _n)] elif _xout = "parameters" then return [seq(_yini[_n+_i], _i = 1 .. nops(_pars))] else return [seq(_yini[_i], _i = 0 .. _n)], [seq(_yini[_n+_i], _i = 1 .. nops(_pars))] end if else return "procname" end if end if; if _xout-_y0[0] = 0. then return [seq(_y0[_i], _i = 0 .. _n)] end if; _reinit := false; if _xin <> "last" then if 0 < 0 and `dsolve/numeric/checkglobals`(0, table( [ ] ), _pars, _n, _yini) then _reinit := true; if _pars <> [] then _par := {seq(rhs(_pars[_i]) = _yini[_n+_i], _i = 1 .. nops(_pars))}; for _i from 0 to _n do _y0[_i] := subs(_par, _yini[_i]) end do; for _i from _n+1 to _n+nops(_pars) do _y0[_i] := _yini[_i] end do else for _i from 0 to _n do _y0[_i] := _yini[_i] end do end if; for _i to _n do _octl[_i+7] := _y0[_i] end do; for _i to nops(_pars) do _octl[2*_n+30+_i] := _y0[_n+_i] end do end if; if _pars <> [] and select(type, {seq(_yini[_n+_i], _i = 1 .. nops(_pars))}, 'undefined') <> {} then error "parameters must be initialized before solution can be computed" end if end if; if not _reinit and _xout-_ctl[2] = 0 then [_ctl[2], seq(_ctl[_i], _i = 8 .. 7+_n)] else if sign(_xout-_ctl[2]) <> sign(_ctl[2]-_y0[0]) or abs(_xout-_y0[0]) < abs(_xout-_ctl[2]) or _reinit then for _i to 34 do _ctl[_i] := _octl[_i] end do end if; _ctl[3] := _xout; if Digits <= evalhf(Digits) then try _errcd := evalhf(`dsolve/numeric/lsode`(_fcn, var(_ctl))) catch: userinfo(2, `dsolve/debug`, print(`Exception in lsode:`, [lastexception])); if searchtext('evalhf', lastexception[2]) <> 0 or searchtext('real', lastexception[2]) <> 0 or searchtext('hardware', lastexception[2]) <> 0 then _errcd := `dsolve/numeric/lsode`(_fcn, _ctl) else error end if end try else _errcd := `dsolve/numeric/lsode`(_fcn, _ctl) end if; if _errcd < 0 then userinfo(2, {dsolve, `dsolve/lsode`}, `Last values returned:`); userinfo(2, {dsolve, `dsolve/lsode`}, ` t =`, _ctl[2]); _i := 8; userinfo(2, {dsolve, `dsolve/lsode`}, ` y =`, _ctl[_i]); for _i from _i+1 to 7+_n do userinfo(2, {dsolve, `dsolve/lsode`}, ` `, _ctl[_i]) end do; if _errcd+1. = 0. then if _ctl[14+trunc(_ctl[1])] <> 0 then error "an excessive amount of work was done, maxstep may be too small" else error "an excessive amount of work (greater than mxstep) was done" end if elif _errcd+2. = 0. then error "too much accuracy was requested for the machine being used" elif _errcd+3. = 0. then error "illegal input was detected" elif _errcd+4. = 0. then error "repeated error test failures on the attempted step" elif _errcd+5. = 0. then error "repeated convergence test failures on the attempted step" elif _errcd+6. = 0. then error "pure relative error control requested for a variable that has vanished" elif _errcd+7. = 0. then error "cannot evaluate the solution past %1, maxfun limit exceeded (see <a href='http://www.maplesoft.com/support/help/search.aspx?term=dsolve,maxfun' target='_new'>?dsolve,maxfun</a> for details)", evalf[8](_ctl[2]) else error "unknown error code returned from lsode %1", trunc(_errcd) end if end if; if _Env_smart_dsolve_numeric = true then if _y0[0] < _xout and procname("right") < _xout then procname("right") := _xout elif _xout < _y0[0] and _xout < procname("left") then procname("left") := _xout end if end if; [_xout, seq(_ctl[_i], _i = 8 .. 7+_n)] end if end proc, (2) = Array(1..3, {(1) = 18446744078356217278, (2) = 18446744078356217454, (3) = 18446744078356217630}), (3) = [x, f[0, 0](x), f[1, 0](x)], (4) = []}); _solnproc := _dat[1]; _pars := map(rhs, _dat[4]); if not type(_xout, 'numeric') then if member(x, ["start", 'start', "method", 'method', "left", 'left', "right", 'right', "leftdata", "rightdata", "enginedata", "eventstop", 'eventstop', "eventclear", 'eventclear', "eventstatus", 'eventstatus', "eventcount", 'eventcount', "laxtol", 'laxtol', "numfun", 'numfun', NULL]) then _res := _solnproc(convert(x, 'string')); if 1 < nops([_res]) then return _res elif type(_res, 'array') then return eval(_res, 1) elif _res <> "procname" then return _res end if elif member(x, ["last", 'last', "initial", 'initial', NULL]) then _res := _solnproc(convert(x, 'string')); if type(_res, 'list') then return _res[3] else return NULL end if elif member(x, ["parameters", 'parameters', "initial_and_parameters", 'initial_and_parameters', NULL]) then _xout := convert(x, 'string'); _res := _solnproc(_xout); if _xout = "parameters" then return [seq(_pars[_i] = _res[_i], _i = 1 .. nops(_pars))] else return [_res[3], seq(_pars[_i] = [_res][2][_i], _i = 1 .. nops(_pars))] end if elif type(_xout, `=`) and member(lhs(_xout), ["initial", 'initial', "parameters", 'parameters', "initial_and_parameters", 'initial_and_parameters', NULL]) then _xout := convert(lhs(x), 'string') = rhs(x); if lhs(_xout) = "initial" then if type(rhs(_xout), 'list') then _res := _solnproc(_xout) else _res := _solnproc("initial" = ["single", 3, rhs(_xout)]) end if elif not type(rhs(_xout), 'list') then error "initial and/or parameter values must be specified in a list" elif lhs(_xout) = "initial_and_parameters" and nops(rhs(_xout)) = nops(_pars)+1 then _res := _solnproc(lhs(_xout) = ["single", 3, op(rhs(_xout))]) else _res := _solnproc(_xout) end if; if lhs(_xout) = "initial" then return _res[3] elif lhs(_xout) = "parameters" then return [seq(_pars[_i] = _res[_i], _i = 1 .. nops(_pars))] else return [_res[3], seq(_pars[_i] = [_res][2][_i], _i = 1 .. nops(_pars))] end if elif type(_xout, `=`) and member(lhs(_xout), ["eventdisable", 'eventdisable', "eventenable", 'eventenable', "eventfired", 'eventfired', "direction", 'direction', NULL]) then return _solnproc(convert(lhs(x), 'string') = rhs(x)) elif _xout = "solnprocedure" then return eval(_solnproc) elif _xout = "sysvars" then return _dat[3] end if; if procname <> unknown then return ('procname')(x) else _ndsol := `tools/gensym`("f[1,0](x)"); eval(FromInert(_Inert_FUNCTION(_Inert_NAME("assign"), _Inert_EXPSEQ(ToInert(_ndsol), _Inert_VERBATIM(pointto(_dat[2][3])))))); return FromInert(_Inert_FUNCTION(ToInert(_ndsol), _Inert_EXPSEQ(ToInert(x)))) end if end if; try _res := _solnproc(_xout); _res[3] catch: error end try end proc

 (4)

ftoplot := unapply(f00(x)+0.45e-1*f10(x)/x^(3/2), x)

proc (x) options operator, arrow; f00(x)+0.45e-1*f10(x)/x^(3/2) end proc

 (5)

``

plot(ftoplot(x), x = xini .. xfin)

  

``


Download Test2ODEs.mw

The approach in that file works, however I have a question regarding the efficiency of my method, since I plan to extend the system to many more ODEs besides just 2 and also extend the range to a larger xfin. In this method, since I define the function to plot in terms of f01 and f02, wich are procedures, does this mean that for each x on the grid for the plot(ftoplot,x=xini..xfin) maple actually computes the solutions f00(x) and f01(x) and then forms the ftoplot combination and plots that specific point? If the default sampling of my interval is, say 1000 points, does it mean that the way I wrote it I will have 1000 invocations of the dsolve procedure, for each x in the sample? I am not sure, it seems to me that is the case. This would imply that instead of advancing the solution at each step maple starts over again from xini. How could I just avoid this behavior and instead have access to the values of ftoplot(x) in the range xini to xfin stored from one invocation of dsolve? 

 

The ideal scenario for me would be to have f[0,0](x) and f[0,1](x) stored as an interpolated function between xini and xfin from the solutions of one invocation of dsolve prior to defining ftoplot. Can this be achieved in principle? How? Remember, i have to use method=lsode and range is not accepted.  

 

Simplify doesn't work for me...

Asked by: h8red 15
April 07 2015
1  4

 

Its pretty obvious that i would like the to get the answer -exp(-n(cos(phi)+sin(phi))) from the 'simplify' command. What am i doing wrong? the problem

Print "Hello" in Green?...

Asked by: mapleatha 160
print
April 07 2015
1  11

How do we print "Hello" in green in a .mw file in Maple 13?

Thanks!

mapleatha

Picking out Christoffel Symbols...

Asked by: shillis 20
April 07 2015
3  2

I would like to pick out specific christoffel symbols once maple calculates them. Right now I am using CoefficientList, but I would like a more direct way, as I am working with someone using Maple 17, and they are put in a different order in the CoefficientList. Here is the bit of my code I would like help with.

Error, invalid input: with expects...

Asked by: bfathi 10
April 07 2015
0  6

hi friends

After this cods i see very error

 > restart;

read(orbit.sav ): whit(plots):
ax := -G*Mz*x/(x^2+y^2)^(3/2);
ay := -G*Mz*y/(x^2+y^2)^(3/2);
i := 'i'; j := i+1;
for k from 0 to 3 do
x := 7*10^6; Vx := 0;
y := 0; Vy := 9000;
dt := evalf(1/2^k);
for i from 0 to 328 do
X[i] := evalf(x); Y[i] := evalf(y);
for n to 40*2^k do
x := evalf((1/2)*ax*dt^2+Vx*dt+x); y := evalf((1/2)*ay*dt^2+Vy*dt+y);
Vx := evalf(ax*dt+Vx); Vy := evalf(ay*dt+Vy)
od;
if i mod 41= 0 then
dX[k, i] := X[i]-XS[j]; dY[k, i] := Y[i]-YS[j]
fi
od;
p[k] := plot([seq([(X[i]-XS[j])*(1/1000), (Y[i]-YS[j])*(1/1000)], i = 0 .. 328)], color = green) end do;
p1 := display({seq(p[k], k = 0 .. 3)}, thickness = 3)
SI := [seq(41*i, i = 0 .. 8)]
p2 := plot({seq([seq([(1/1000)*dX[k, i], (1/1000)*dY[k, i]], k = 0 .. 1), [0, 0]], i = SI)}, color = black)
display({p1, p2}, scaling = constrained, labels = ['dx', 'dy'])
display({p1, p2}, view = [-.1 .. .5, -.4 .. .2], scaling = constrained, labels = ['dx', 'dy'])

can you help me Please?

Thank you

 

 

 

ColumnGraph in 3D...

Asked by: Alger 524
3d column plotting
April 07 2015
1  7

Dear all,

Is there a Maple command to get ColumnGraph in 3D as bar3 in Matlab.

In 2D, I use:

Statistics[ColumnGraph](A);  # equivalent to bar(A) in Matlab

Thanks

how do i solve this problem?...

Asked by: hooshmand 5
fsolve
April 07 2015
1  1

Hi,

we want to know what is the meaning of this statement?

We expect to have the following statement . But unfortunately we don,t get it

> with(difforms);
> sol := fsolve({diff(S, x) = 0, diff(S, y) = 0}, {x, y});


I do not take values above code.

Can't figure this out...

Asked by: Rouben Rostamian 8531
April 06 2015
1  7

Consider

> z := Int(f(t-s), s=0..1);

Int(f(t-s), s = 0 .. 1)

The result is clearly a function of t (and definitely not a function of s, which is a dummy variable).  Let's define

> F := unapply(z, t);

proc (t) options operator, arrow; Int(f(t-s), s = 0 .. 1) end proc

Then

> F(q);

Int(f(q-s), s = 0 .. 1)

as expected.  However:

> F(s);

Int(f(0), s = 0 .. 1)

which is not correct, as this confuses the argument s of F with the dummy integration variable s.  How would you salvage the situation?

 

Problem with integration....

Asked by: schzan 130
April 06 2015
0  11

I am performing the int operation on an experssion but it does not consider constants to be constants. For example, maple is considering this constant namely, f'(-1) or eval(diff(f(x), x), x=-1) as a function and henceforth does not give the desired output. Please tell me how to make eval(diff(f(x), x), x=-1) as a constant function or to say, how to convert any given function in to a constant function.

Thank you for your time and help!

Error, invalid left hand side in assignment...

Asked by: wenny 25
April 06 2015
1  2

Here is my code

ma[1] := [-885.880598855072776, [bh = 0., g0h = 0., g1h = 0.825946224722250e-4]]

ma[2] := [-877.957885609114328, [bh = 347.116836805625, g0h = 0., g1h = 0.164861392564e-3]]

ma[3] := [-863.445203144006655, [bh = 0., g0h = 0., g1h = 0.787090377403668e-4]]

avb := 90.5276611129000; avg0 := 0; avg1 := 0.92225359e-4;

for j from 1 to 3 do

assign(ma[j][2]);

A[j] := [(avb-bh)^2, (avg0-g0h)^2, (avg1-g1h)^2] 

end do;

 

Result is that

A[1] := [4763.19965962732, 0., 1.13103562500664*10^(-10)]

Error, invalid left hand side in assignment

 it only shows A[1], but no A[2] and A[3]

Please help!

Convolution animation...

Asked by: nMaple 85
animation
April 06 2015
1  21

Hello,

Is it possible to create animation of convolution of two functions?

For example f(t)=u(t)-u(t-2) and g(t) = tu(t)-(t-4)u(t-4), where u(t) is a step function.

I would like to generate animation for this convolution.

Any help would be appreciated.

 

Thanks.

My math discovery done with Maple

by: Markiyan Hirnyk 7819 Maple
cases
April 06 2015
4  2

 

     Maple is seriously used in my article Approximation of subharmonic functions in the half-plane by the logarithm of the modulus of an analytic function. Math. Notes 78, No 4, 447-455 in two places. The purpose of this post is to present these applications.                                                                                                 First, I needed to prove the elementary inequality (related to the properties of the minimal harmonic majorant of the function 1/Im z in a certain strip)                                                                                                    2R+sqrt(R)-R(R+sqrt(R))y - 1/y ≤  1/4                                                                                                  for    y ≥ 1/(R+sqrt(R)) and  y ≤ 1/R, the parameter R is greater than or equal to 1.   The artless attemt                                                                          
restart; `assuming`([maximize(2*R+sqrt(R)-R*(R+sqrt(R))*y-1/y, y = 1/(R+sqrt(R)) .. 1/R)], [R >= 1])

maximize(2*R+R^(1/2)-R*(R+R^(1/2))*y-1/y, y = 1/(R+R^(1/2)) .. 1/R)

 (1)

fails. The second (and successful) try consists in the use of optimizers:

F := proc (R) options operator, arrow; evalf(maximize(2*R+sqrt(R)-R*(R+sqrt(R))*y-1/y, y = 1/(R+sqrt(R)) .. 1/R)) end proc:

F(1)

.171572876

 (2)

 

Optimization:-Minimize('F(R)', {R >= 1})

[.171572875253809986, [R = HFloat(1.0)]]

 (3)

To be sure ,
DirectSearch:-Search(proc (R) options operator, arrow; F(R) end proc, {R >= 1})
;

[.171572875745665, Vector(1, {(1) = 1.0000000195752754}, datatype = float[8]), 11]

 (4)

Because 0.17
"158 < 0.25, the inequality is proved. "
Now we establish this  by the use of the derivative. 

solve(diff(2*R+sqrt(R)-R*(R+sqrt(R))*y-1/y, y) = 0, y, explicit)

1/(R^(3/2)+R^2)^(1/2), -1/(R^(3/2)+R^2)^(1/2)

 (5)

maximize(1/sqrt(R^(3/2)+R^2)-1/(R+sqrt(R)), R = 1 .. infinity, location)

(1/2)*2^(1/2)-1/2, {[{R = 1}, (1/2)*2^(1/2)-1/2]}

 (6)

minimize(eval(2*R+sqrt(R)-R*(R+sqrt(R))*y-1/y, y = 1/sqrt(R^(3/2)+R^2)), R = 1 .. infinity, location)

3-2*2^(1/2), {[{R = 1}, 3-2*2^(1/2)]}

 (7)

evalf(3-2*sqrt(2))

.171572876

 (8)

The second use of Maple was the calculation of the asymptotics of the following integral (This is the double integral of the Laplacian of 1/Im z over the domain {z: |z-iR/2| < R/2} \ {z: |z| ≤ 1}.). That place is the key point of the proof. Its direct calculation in the polar coordinates fails.

`assuming`([(int(int(2/(r^2*sin(phi)^3), r = 1 .. R*sin(phi)), phi = arcsin(1/R) .. Pi-arcsin(1/R)))/(2*Pi)], [R >= 1])

(1/2)*(int(int(2/(r^2*sin(phi)^3), r = 1 .. R*sin(phi)), phi = arcsin(1/R) .. Pi-arcsin(1/R)))/Pi

 (9)

In order to overcome the difficulty, we find the inner integral

`assuming`([(int(2/(r^2*sin(phi)^3), r = 1 .. R*sin(phi)))/(2*Pi)], [R*sin(phi) >= 1])

(R*sin(phi)-1)/(sin(phi)^4*R*Pi)

 (10)

and then we find the outer integral. Because
`assuming`([int((R*sin(phi)-1)/(sin(phi)^4*R*Pi), phi = arcsin(1/R) .. Pi-arcsin(1/R))], [R >= 1])

int((R*sin(phi)-1)/(sin(phi)^4*R*Pi), phi = arcsin(1/R) .. Pi-arcsin(1/R))

 (11)

is not successful, we find the indefinite integral  

J := int((R*sin(phi)-1)/(sin(phi)^4*R*Pi), phi)

-(1/2)*cos(phi)/(Pi*sin(phi)^2)+(1/2)*ln(csc(phi)-cot(phi))/Pi+(1/3)*cos(phi)/(R*Pi*sin(phi)^3)+(2/3)*cos(phi)/(R*Pi*sin(phi))

 (12)

We verify that  the domain of the antiderivative includes the range of the integration.
plot(-cos(phi)/sin(phi)^2+ln(csc(phi)-cot(phi)), phi = 0 .. Pi)

  

plot((2/3)*cos(phi)/sin(phi)^3+(4/3)*cos(phi)/sin(phi), phi = 0 .. Pi)

  

    That's all right. By the Newton-Leibnitz formula,

``
eval(J, phi = Pi-arcsin(1/R))-(eval(J, phi = arcsin(1/R)));

(1/3)*(1-1/R^2)^(1/2)*R^2/Pi+(1/2)*ln((1-1/R^2)^(1/2)*R+R)/Pi-(4/3)*(1-1/R^2)^(1/2)/Pi-(1/2)*ln(R-(1-1/R^2)^(1/2)*R)/Pi

 (13)

Finally, the*asymptotics*is found by

asympt(eval(J, phi = Pi-arcsin(1/R))-(eval(J, phi = arcsin(1/R))), R, 3)

(1/3)*R^2/Pi-(3/2)/Pi+(1/2)*(ln(2)+ln(R))/Pi-(1/2)*(-ln(2)-ln(R))/Pi+O(1/R^2)

 (14)

      It should be noted that a somewhat different expression is written in the article. My inaccuracy, as far as I remember it, consisted in the integration over the whole disk {z: |z-iR/2| < R/2} instead of {z: |z-iR/2| < R/2} \ {z: |z| ≤ 1}. Because only the form of the asymptotics const*R^2 + remainder is used in the article, the exact value of this non-zero constant is of no importance.

       It would be nice if somebody else presents similar examples here or elsewhere.

 

Download Discovery_with_Maple.mw

First 1250 1251 1252 1253 1254 1255 1256 Last Page 1252 of 2224