:-)
even if that sounds funny, it is a big mess (and I am quite sure it violates users rights in the EC)
practically I have several Excel sheets and I *have* to test them in future time (for correct behaviour w.r.t. calendars)
otherwise said: not only Bill Gates seems to believe each computer his is own ...
I see it it is a problem if you need a new PC.
But may be you check some info pages before, say Wikipedia - where i found the German site quite instructive, http://de.wikipedia.org/wiki/Microsoft_Windows_Vista. And it sounds horrible: you pay and get no rights for that ...
And US Departments forbid to upgrade http://www.dot.gov/ost/m60/morat001.pdf
Time to think about *nix and whether & how old Office software can be used.
No series expansion needed I think, it should be all in the context
of finite series, the rest are formal embeddings.
Primarily it is not a Maple problem and IIRC it goes like this:
f smooth in 0 implies there is a w in the ring s.th. k[x,y] = k[w,f]
is a local isomorphism and it induces ring/(f) = k[w] in 0 (this is what
you call the embedding, but it is achieved that way).
Essentially this follows from the implicite function theorem.
Geometrically it means: if the curve is smooth you have a (local) tangent,
a complementary vector spans the space. Which are your new coordinates.
One has to be a bit carefully since you are in characteristic p > 0.
For the final step - Laurent series - be careful, not sure whether Maple
needs 1/n! which is not given in your case. What you mean is the completed
local ring (a limit construction) to which the stuff always embeds. Then
the Laurent series are the quotient field (as your local ring is integer).
I am not aware of a directly applicable package, may be you dig for
'abstract algebra' and Maple, may be http://mihailovs.com/Alec/ .
Or other CAS software dedicated to commutative algebra / algebraic geometry.
But if you write down the Math then it should translatable to Maple later,
i think it allows finite fields (do not know whether abstract or only for
explicit cases).
Roughly prices are log-normal (not normal, this would be Bachelier I think),
may be you google for lectures using 'binomial option pricing' (or similar)
[but be warned: social reality has no actual reason to follow ideal models].
Doing a search at this board for 'brownian' you will find some nice stuff
for that on Alec's blog for simulations.
If you use Compile it will become even faster (as generating the compiled version
will not count - and that's what would be done in VB as well).
BTW: VB should give you ~ 15 places (double precision), i think you have to use
type conversion before output, something like (cdbl)total/(cdbl)1000000.
I hate mw sheets ... you fill your pq with integer values between 1 and 70 (it is an array with 2 columns),
that may be done certainly faster, but it does not kill you:
Even if you start with those integers I think through H1,H2,H3 you will pass to rationals (have not checked),
so it may be worth to consider working with floats (i.e. with evalhf).
Your matrix Theta is not initialized, so it is filled with 0 and Vpq[s,i] := Rpq[i,s]*exp(I*pq[s,2]*Theta[i,1])
will simply be Rqp, but transposed.
If this is a simple coding error and it should have a reasonable value at least here you will fall out of
the rationals.
And as it seems you finally want the vector 'moments' but you do not return it (but do not understand,
why you need 'Typesetting:-delayDotProduct' ??)
So I would suggest not to work with integers and symbolic exactness, but to work with pure floating
point computations, which is much faster and may even be used in compiled form. It may be that
the redim could make some trouble for it, but you actually do not need it (just define the matrices
in the form you want them ... sounds as if you did a lot of Excel coding ...)
Hope I got it right ... I have not looked at in detail is your last double loop, it is of magnitude
size * dim^2 ~ 1/4*(pmax+2)^2 * 200 ^2 ~ 72^2/4 * 200^2 = 51.840.000 of body calls and you are
silently killed by a combinatorical explosions I guess - how fast ever the numerics will be: this
may be too much for anybodys patients.
As I stumbled over a similar problem these days here is a suggestion.
The recursion can be made linear by z(t) = y(t+1)/y(t):
# make it linear by substitution
z(t+1)=1+a*b-((a*b)/z(t)) :
subs( z(t) = y(t+1)/y(t), %):
subs( z(t+1) = y(t+2)/y(t+1), %):
isolate(%, y(t+2)):
# bring it in nice shape (not needed)
expand(%): collect(%, y(t+1));
y(t + 2) = (1 + a b) y(t + 1) - a b y(t)
Now solve the new 2-term linear recursion and assign the solution
rsolve({%, y(0)=alpha, y(1)=beta}, y(k)):
Y:=unapply(%,k):
'Y(k)': '%'=%;
k
(alpha - beta) (a b) -a b alpha + beta
Y(k) = - --------------------- - -----------------
a b - 1 a b - 1
From that we get a general solution for z
z(t) = 'eval(y(t+1)/y(t), y=Y)':
simplify(%): rhs(%):
Z:=unapply(%,t):
'Z(t)': '%'=%;
t t
(a b) a b alpha - (a b) a b beta - a b alpha + beta
Z(t) = - -----------------------------------------------------
t t
-(a b) alpha + (a b) beta + a b alpha - beta
Now we have to recover z from Z for which one has to determine the indeterminates.
As we have 2 of them we want 2 equations, let us take t=0 and t=1:
'Z(0)': '%'= simplify(%);
'Z(1)': '%'= collect(collect(simplify(%),a),b);
beta
Z(0) = -----
alpha
(alpha - beta) a b
Z(1) = 1 - ------------------
beta
Solving the first equation for beta = z0*alpha; and inserting that into the
second equation we get no condition on alpha as it cancels out.
So set alpha to 1 and a solution is found with beta = z(0):
z(k) = eval(Z(k), [alpha=1, beta=z(0)]);
k k
(a b) a b - (a b) a b z(0) - a b + z(0)
z(k) = - -----------------------------------------
k k
-(a b) + (a b) z(0) + a b - z(0)
you will want to look at the circle around 0 - 3/2*I of radius sqrt(13)/2 and watch
the angles at arctan(3/2) and Pi-arccot(2/3), that should be your task (I think ...):
K0 := -3/2*I+1/2*13^(1/2)*exp(t*I);
plot( [Re(K0), Im(K0), t=arctan(3/2) .. Pi-arccot(2/3)], scaling=constrained, color=blue):
plot( [Re(K0), Im(K0), t=-Pi .. arctan(3/2)], scaling=constrained, color=red):
plot( [Re(K0), Im(K0), t=Pi-arccot(2/3) .. Pi], scaling=constrained, color=magenta):
plots[display](%,%%,%%%);
i think maple cleans up blanks, so it removes any you add
f:=x/(1/2-x^2)^(1/3);
F:=-3/4*(1/2-x^2)^(2/3);
'diff(F,x) = f'; is(%);
x
f := -------------
2 1/3
(1/2 - x )
2 2/3
3 (1/2 - x )
F := - ---------------
4
d
-- F = f
dx
true
'[Int(f, x = 0 .. 1/2*2^(1/2)), Int(f, x = 1/2*2^(1/2) .. 1)]';
value(%);
evalf(%);
%[1] + %[2];
1/2
2
----
2 1
/ /
| |
| f dx, | f dx
| |
/ /
0 1/2
2
----
2
1/3 1/3
3 2 3 2 (1/3) 1/2
[------, ------ - 3/16 I 2 3 ]
8 16
[0.4724703938, 0.2362351969 - 0.4091713637 I]
0.7087055907 - 0.4091713637 I
['limit(F,x=1/2*2^(1/2), left) - eval(F,x=0)', 'eval(F,x=1) - limit(F,x=1/2*2^(1/2), right)'];
%;
evalf(%);
%[1] + %[2];
/ lim F\ - F| , F| - / lim F\
| / 1/2\ | |x = 0 |x = 1 | / 1/2\ |
| |2 | | | |2 | |
|x -> |----|- | |x -> |----|+ |
\ \ 2 / / \ \ 2 / /
1/3 2/3 1/3
3 2 3 (-1) 2
[------, - --------------]
8 8
[0.4724703938, 0.2362351969 - 0.4091713635 I]
0.7087055907 - 0.4091713635 I
Int(f,x=0..1);
value(%);
evalf(%);
1
/
| x
| ------------- dx
| 2 1/3
/ (1/2 - x )
0
1/3
9 2 (1/3) 1/2
------ - 3/16 I 2 3
16
0.7087055906 - 0.4091713637 I
# good practice: restart your system to clean up
restart;
# your mag is usually called the absolute value
sqrt(z*conjugate(z))- abs(z); simplify(%);
_ 1/2
(z z) - | z |
0
Zeta(0,n+1,1-2*Pi*v*I)/Zeta(n+1):
indets(%,atomic): indets(%); # shows you all indeterminates
{v, n}
# define it as function
f := (n,v) -> Zeta(0,n+1,1-2*Pi*v*I)/Zeta(n+1);
Zeta(0, n + 1, 1 - 2 I Pi v)
f := (n, v) -> ----------------------------
Zeta(n + 1)
# provide some n
nTst:= 5;
nTst := 5
# set up and solve the equation with 14 Digits
Digits:=14;
abs( f(nTst,v) ) = 0.5;
theSol:=fsolve(%,v);
Digits := 14
| Psi(5, 1 - 2 I Pi v) |
63/8 ------------------------ = 0.5
6
Pi
theSol := 0.019705381874366 I
# test the solution
f(nTst, theSol): evalf(%);
0.49999999999993
Not sure what you expect as an answer ... I used the following to re-write it
subs(n1=k,Ans);
subs(Q(0)=Q0,%);
ans:=subs(n=m+1,%);
and then used assume(k::nonnegint); additionally(k LE m); getassumptions( [k,m]);
(read LE as less then, Will should really give some option pure text with fixed font)
Now simplify(ans); combine(%); gives an expression, where you only need the Sum
and kicking off summation constants I achieve at the following:
Sum((2^(2*m)*4^(-k)+3*4^(m-k)*k+2*4^(m-k)*k^2) /
GAMMA((k*x+1/2*(-x)^(1/2)+2*x)/x)/GAMMA((k*x-1/2*(-x)^(1/2)+2*x)/x)*k,k = 0 .. m);
The numerator is 4^(m-k)*(2*k^2+3*k+1)*k and the denomerator is a product of 2
Gamma functions Gamma(k-u)*Gamma(k+u), 1/2*(-x)^(1/2)/x+2 = u.
That is something involving sin or sinc or may be written through Pochhammer I
think and now I ask Maple again for an answer:
4^(m-k)*(2*k^2+3*k+1)*k / GAMMA(k-u)/GAMMA(k+u); Sum(%,k=0..m); value(%);
(m - k) 2
4 (2 k + 3 k + 1) k
---------------------------
GAMMA(k - u) GAMMA(u + k)
m
----- (m - k) 2
\ 4 (2 k + 3 k + 1) k
) ---------------------------
/ GAMMA(k - u) GAMMA(u + k)
-----
k = 0
m
-3/2 hypergeom([5/2, 3], [3/2, 1 - u, u + 1], 1/4) 4
sin(Pi (u + 1))/(u Pi) - 1/2 (m + 3/2) (m + 2) (m + 1)
hypergeom([1, m + 3, m + 5/2],
[m + 1, u + m + 1, m + 1 - u, m + 3/2], 1/4)/(
GAMMA(m + 1 - u) GAMMA(u + m + 1))
As Maple has some bugs on hypergeometrics take it with care and check against
the case where m -> infinity:
4^(m-k)*(2*k^2+3*k+1)*k / GAMMA(k-u)/GAMMA(k+u): Sum(%,k=0..infinity); value(%);
infinity
----- (m - k) 2
\ 4 (2 k + 3 k + 1) k
) ---------------------------
/ GAMMA(k - u) GAMMA(u + k)
-----
k = 0
m
-3/2 hypergeom([5/2, 3], [3/2, 1 - u, u + 1], 1/4) 4
sin(Pi (u + 1))/(u Pi)
I hope I made not too much errors and you can build your answer from that:
the answer is essential: you have the finite summation part of some 2F3 * sinc
I see no way to get it in 'simple' form. If it is a practical problem, the you may know even more about possible solutions looking at data and using numerics - or by analyzing the original problem. But a general answer may often be lengthy (like the above)
Possibly I do not quite understand what you say, it is an analytic solution
(up to degree 4 there are formulae in terms of roots and i think maple just
uses that).
sols[1]: indets(%, atomic): indets(%);
{L, x1}
This shows, that the expression (the first of 4 solutions) depends on x1 and
L and if you want to make it a function use the following command:
f1 := unapply( sols[1], x1,L):
Then you can use it as a function (it is about 12 pages long, 60867 characters
80 per line and 60 lines per page ~ 12).
If you want it purely numerical, then you need numerical values for x1 and L.
The reduction using the range you know seems not to work, may be, that your x1
and L are too arbitrary - but the expression is too long to expect to much.
What is the reason for your question, what is the actual problem?
