Maple Questions and Posts

These are Posts and Questions associated with the product, Maple

Hello,

I have code which contains a data structure. This data structure was defined with Types and a Table.

I would like to know which is the best manner to do to build a data structure.

Is it better to use a Table or a MultiSet ?

May you precise me the pros and the cons of the 2 solutions?

Or indicate me if you have other ways to do.

Thanks a lot for your help.

I want to sum a series by getting multiple values from a list.  Based on the list L I was trying

add(A, [m2, m3] in L)   This seem to have dissapeared from the document but the error message remains.

restart

``

NULL

NULL

A := factorial(2*m2+3*m3)*c0^(1+m2+2*m3)*c2^m2*c3^m3/(factorial(1+m2+2*m3)*factorial(m2)*factorial(m3)*c1^(2*m2+3*m3+1))

factorial(2*m2+3*m3)*c0^(1+m2+2*m3)*c2^m2*c3^m3/(factorial(1+m2+2*m3)*factorial(m2)*factorial(m3)*c1^(2*m2+3*m3+1))

(1)

"P:=3, N:=4"

4

(2)

``

I found this code by Kitonum 16740  https://www.mapleprimes.com/questions/138955-How-Do-I-Create-A-Variable-Number-Of-Nested-Loops-

S := {seq(0 .. N)}; l := P-1; f := proc (T) options operator, arrow; [seq(seq([op(T[j]), S[i]], i = 1 .. nops(S)), j = 1 .. nops(T))] end proc; L := (f@@l)([[]])

[[0, 0], [0, 1], [0, 2], [0, 3], [0, 4], [1, 0], [1, 1], [1, 2], [1, 3], [1, 4], [2, 0], [2, 1], [2, 2], [2, 3], [2, 4], [3, 0], [3, 1], [3, 2], [3, 3], [3, 4], [4, 0], [4, 1], [4, 2], [4, 3], [4, 4]]

(3)

nops(L)

25

(4)

NULLNULL

``

nops(L)

25

(5)

Instead of using add(add(A, m2 = 0 .. 4)*m3 = 0 .. 4)I would like to map m2 and m3 to the values in list L and sum them. This should give a method for having a variable number of loops

I tried variants off without success

 

add(A, `in`([m2, m3], L))

Error, illegal use of an object as a name

 

 

 

 

 

for the next higher order would need to map m2,m3, m4 to the generated list.

 

 

 

B := factorial(2*m2+3*m3+4*m4)*c0^(1+m2+2*m3+3*m4)*c2^m2*c3^m3*c4^m4/(factorial(1+m2+2*m3+3*m4)*factorial(m2)*factorial(m3)*factorial(m4)*c1^(2*m2+3*m3+4*m4+1))

factorial(2*m2+3*m3+4*m4)*c0^(1+m2+2*m3+3*m4)*c2^m2*c3^m3*c4^m4/(factorial(1+m2+2*m3+3*m4)*factorial(m2)*factorial(m3)*factorial(m4)*c1^(2*m2+3*m3+4*m4+1))

(6)

P := 4; S := {seq(0 .. N)}; l := P-1; f := proc (T) options operator, arrow; [seq(seq([op(T[j]), S[i]], i = 1 .. nops(S)), j = 1 .. nops(T))] end proc; L := (f@@l)([[]])

[[0, 0, 0], [0, 0, 1], [0, 0, 2], [0, 0, 3], [0, 0, 4], [0, 1, 0], [0, 1, 1], [0, 1, 2], [0, 1, 3], [0, 1, 4], [0, 2, 0], [0, 2, 1], [0, 2, 2], [0, 2, 3], [0, 2, 4], [0, 3, 0], [0, 3, 1], [0, 3, 2], [0, 3, 3], [0, 3, 4], [0, 4, 0], [0, 4, 1], [0, 4, 2], [0, 4, 3], [0, 4, 4], [1, 0, 0], [1, 0, 1], [1, 0, 2], [1, 0, 3], [1, 0, 4], [1, 1, 0], [1, 1, 1], [1, 1, 2], [1, 1, 3], [1, 1, 4], [1, 2, 0], [1, 2, 1], [1, 2, 2], [1, 2, 3], [1, 2, 4], [1, 3, 0], [1, 3, 1], [1, 3, 2], [1, 3, 3], [1, 3, 4], [1, 4, 0], [1, 4, 1], [1, 4, 2], [1, 4, 3], [1, 4, 4], [2, 0, 0], [2, 0, 1], [2, 0, 2], [2, 0, 3], [2, 0, 4], [2, 1, 0], [2, 1, 1], [2, 1, 2], [2, 1, 3], [2, 1, 4], [2, 2, 0], [2, 2, 1], [2, 2, 2], [2, 2, 3], [2, 2, 4], [2, 3, 0], [2, 3, 1], [2, 3, 2], [2, 3, 3], [2, 3, 4], [2, 4, 0], [2, 4, 1], [2, 4, 2], [2, 4, 3], [2, 4, 4], [3, 0, 0], [3, 0, 1], [3, 0, 2], [3, 0, 3], [3, 0, 4], [3, 1, 0], [3, 1, 1], [3, 1, 2], [3, 1, 3], [3, 1, 4], [3, 2, 0], [3, 2, 1], [3, 2, 2], [3, 2, 3], [3, 2, 4], [3, 3, 0], [3, 3, 1], [3, 3, 2], [3, 3, 3], [3, 3, 4], [3, 4, 0], [3, 4, 1], [3, 4, 2], [3, 4, 3], [3, 4, 4], [4, 0, 0], [4, 0, 1], [4, 0, 2], [4, 0, 3], [4, 0, 4], [4, 1, 0], [4, 1, 1], [4, 1, 2], [4, 1, 3], [4, 1, 4], [4, 2, 0], [4, 2, 1], [4, 2, 2], [4, 2, 3], [4, 2, 4], [4, 3, 0], [4, 3, 1], [4, 3, 2], [4, 3, 3], [4, 3, 4], [4, 4, 0], [4, 4, 1], [4, 4, 2], [4, 4, 3], [4, 4, 4]]

(7)

nops(L)

125

(8)

NULL


 

Download Sum_by_mapping_to_a_list_question.mw

sum2N := proc (N::integer) local summation, i, k; summation := 0; k := 0;

for i from 0 to N do if i = 1 then k := k+1

end if; summation := (i+k)^2 end do

end proc;
sum2N(10);

the result shows 121 but its not correct i think the right answer show be 2985
 

Bonjour,

Comment transformer cette liste de résultat: 

a := {x = -12, y = -23}, {x = -34, y = 87}, {x = 18, y = 98}, {x = 27, y = 93}, {x = 45, y = 68}

En une liste simple telle que L: = [-12,-23,-34,87,18,98,27,93,45,68 ]

Merci 


HELLO,

How to transform this result list : 

a := {x = -12, y = -23}, {x = -34, y = 87}, {x = 18, y = 98}, {x = 27, y = 93}, {x = 45, y = 68}

In a simple list such that :  L: = [-12,-23,-34,87,18,98,27,93,45,68 ]

Thank You

I have two summations that Maple converts to hypergeometric outputs. Is there a way to stop that?

I dont know how to get Maple to expand the Hyper geometric output for a given walue of N

if the value of N is defined first the output is as I would expect.

The real problem, the sumations are nested and when I get to four nested sums Maple sits there for hours trying to evaluate the general case. 

I have only included the first 2 summations.

Any insights on this would be appreciated.
 

restart

NULL

NULL

NULL

``

``

a := -(sum(C2^(m-1)*C0^m*factorial(2*m)/(C1^(2*m-1)*factorial(m+1)*factorial(m)), m = 0 .. N))

-2*C1/(C2*(1+(-(4*C0*C2-C1^2)/C1^2)^(1/2)))+C2^N*C0^(N+1)*GAMMA(2*N+3)*(N+2)*hypergeom([1, 3/2+N], [N+3], 4*C0*C2/C1^2)/(C1^(2*N+1)*GAMMA(N+3)^2)

(1)

N := 4

4

(2)

a

-2*C1/(C2*(1+(-(4*C0*C2-C1^2)/C1^2)^(1/2)))+42*C2^4*C0^5*hypergeom([1, 11/2], [7], 4*C0*C2/C1^2)/C1^9

(3)

b := -(sum(C2^(m-1)*C0^m*factorial(2*m)/(C1^(2*m-1)*factorial(m+1)*factorial(m)), m = 0 .. N))

-C1/C2-C0/C1-2*C2*C0^2/C1^3-5*C2^2*C0^3/C1^5-14*C2^3*C0^4/C1^7

(4)

NULL

NULL``

unassign('N')

c := sum(sum((-1)^(m[3]+1)*factorial(2*m[2]+3*m[3])*C0^(1+m[2]+2*m[3])*C2^m[2]*C3^m[3]/(factorial(1+m[2]+2*m[3])*factorial(m[2])*factorial(m[3])*C1^(1+2*m[2]+3*m[3])), m[3] = 0 .. N), m[2] = 0 .. N)

sum(C2^m[2]*(-factorial(2*m[2])*C0^(1+m[2])*hypergeom([(2/3)*m[2]+1, (2/3)*m[2]+2/3, 1/3+(2/3)*m[2]], [1+(1/2)*m[2], 3/2+(1/2)*m[2]], -(27/4)*C0^2*C3/C1^3)/(factorial(1+m[2])*C1^(1+2*m[2]))-factorial(2*m[2]+3*N+3)*C0^(3+m[2]+2*N)*C3^(N+1)*hypergeom([1, (2/3)*m[2]+2+N, 5/3+(2/3)*m[2]+N, 4/3+(2/3)*m[2]+N], [N+2, 2+(1/2)*m[2]+N, 5/2+(1/2)*m[2]+N], -(27/4)*C0^2*C3/C1^3)*(-1)^N/(factorial(3+m[2]+2*N)*factorial(N+1)*C1^(4+2*m[2]+3*N)))/factorial(m[2]), m[2] = 0 .. N)

(5)

NULL

NULLN := 4

4

(6)

c

-C0*hypergeom([1/3, 2/3], [3/2], -(27/4)*C0^2*C3/C1^3)/C1-273*C0^11*C3^5*hypergeom([1, 16/3, 17/3], [6, 13/2], -(27/4)*C0^2*C3/C1^3)/C1^16+C2*(-C0^2*hypergeom([1, 4/3, 5/3], [3/2, 2], -(27/4)*C0^2*C3/C1^3)/C1^3-6188*C0^12*C3^5*hypergeom([1, 19/3, 20/3], [13/2, 7], -(27/4)*C0^2*C3/C1^3)/C1^18)+(1/2)*C2^2*(-4*C0^3*hypergeom([5/3, 7/3], [5/2], -(27/4)*C0^2*C3/C1^3)/C1^5-162792*C0^13*C3^5*hypergeom([1, 20/3, 22/3], [6, 15/2], -(27/4)*C0^2*C3/C1^3)/C1^20)+(1/6)*C2^3*(-30*C0^4*hypergeom([7/3, 8/3], [5/2], -(27/4)*C0^2*C3/C1^3)/C1^7-4883760*C0^14*C3^5*hypergeom([1, 22/3, 23/3], [6, 15/2], -(27/4)*C0^2*C3/C1^3)/C1^22)+(1/24)*C2^4*(-336*C0^5*hypergeom([10/3, 11/3], [7/2], -(27/4)*C0^2*C3/C1^3)/C1^9-164745504*C0^15*C3^5*hypergeom([1, 25/3, 26/3], [6, 17/2], -(27/4)*C0^2*C3/C1^3)/C1^24)

(7)

NULL

e := sum(sum((-1)^(m[3]+1)*factorial(2*m[2]+3*m[3])*C0^(1+m[2]+2*m[3])*C2^m[2]*C3^m[3]/(factorial(1+m[2]+2*m[3])*factorial(m[2])*factorial(m[3])*C1^(1+2*m[2]+3*m[3])), m[3] = 0 .. N), m[2] = 0 .. N)

-C0/C1-5*C2^3*C0^4/C1^7-C2*C0^2/C1^3-2*C2^2*C0^3/C1^5-3*C0^5*C3^2/C1^7+12*C0^7*C3^3/C1^10-55*C0^9*C3^4/C1^13-14*C0^5*C2^4/C1^9+C0^3*C3/C1^4+330*C0^7*C2^4*C3/C1^12-5005*C0^9*C2^4*C3^2/C1^15+61880*C0^11*C2^4*C3^3/C1^18-678300*C0^13*C2^4*C3^4/C1^21+5*C0^4*C2*C3/C1^6-28*C0^6*C2*C3^2/C1^9+165*C0^8*C2*C3^3/C1^12-1001*C0^10*C2*C3^4/C1^15+21*C0^5*C2^2*C3/C1^8-180*C0^7*C2^2*C3^2/C1^11+1430*C0^9*C2^2*C3^3/C1^14-10920*C0^11*C2^2*C3^4/C1^17+84*C0^6*C2^3*C3/C1^10-990*C0^8*C2^3*C3^2/C1^13+10010*C0^10*C2^3*C3^3/C1^16-92820*C0^12*C2^3*C3^4/C1^19

(8)

NULL


 

Download why_hypergeom.mw

I am a student who is learning to use Maple to perform numerical analysis of phase field. There has been no solution to the problem of high anisotropy. The reference code is as follows

W :=piecewise(beta>=Pi/2+thetaM and beta<=Pi-thetaM, 1+eps*cos(4*beta), 
       beta>=Pi/2-thetaM and beta<=Pi/2+thetaM, (1+esp*cos(4*thetaM))/cos(thetaM)*cos(4*beta));

 

However, there is no problem at all when simply performing the following settings.

W := 1+ eps*cos(4*beta);

Please tell me how to write such a conditional maple statement. I have referenced many examples on the official website, but none of them succeeded. This includes writing the code in the following form.

W :=betta->piecewise(beta>=Pi/2+thetaM and beta<=Pi-thetaM, 1+eps*cos(4*beta), 
       beta>=Pi/2-thetaM and beta<=Pi/2+thetaM, (1+esp*cos(4*thetaM))/cos(thetaM)*cos(4*beta));

Remarks, beta is an angle, set as follows.

subs(beta=t_atan2(diff(pho(x,y),y), diff(pho(x,y),x))

 


 

We give a line (D) and a point A located at a distance AH=h from D. A constant angle of magnitude alpha pivots to its apex A and we call B and C the points where its sides cut the line D. Let O be the center of the circle circumscribed to the triangle ABC.
Demonstrate that the B and C tangents to the O circle keep a fixed direction. 
Here is my code which don't work for slopes are not equal.

restart; with(plots): with(geometry):unprotect(D):
_EnvHorizontalName := 'x':_EnvVerticalName := 'y':
line(D, y = (1/2)*x-1); point(A, 5, 5); PerpendicularLine(lp, A, D); h := distance(A, D); intersection(H, D, lp);
alpha := (1/16)*Pi;
rotation(lp1, lp, (1/6)*Pi, 'clockwise', A); rotation(lp2, lp1, (1/6)*Pi-alpha, 'clockwise', A); FindAngle(lp1, lp2); evalf(%);
intersection(B, D, lp1); intersection(C, D, lp2);
triangle(T, [A, B, C]);
circumcircle(Elc, T, 'centername' = OO);
TangentLine(tgB, B, Elc); TangentLine(tgC, C, Elc);
evalf(slope(tgB)); evalf(slope(tgC));
dr := draw([D(color = blue), lp(color = red), Elc(color = green), A, B, C, T(color = black), H, tgB, tgC], printtext = true);

display([dr], axes = none, scaling = constrained);
Fig := proc (k) local dr, Elc, B, C, lp1, lp2; global D, A, lp, H, alpha; geometry:-rotation(lp1, lp, (1/6)*Pi+k, 'clockwise', A); geometry:-rotation(lp2, lp1, (1/6)*Pi-alpha+k, 'clockwise', A); geometry:-intersection(B, D, lp1); geometry:-intersection(C, D, lp2); geometry:-triangle(T, [A, B, C]); geometry:-circumcircle(Elc, T, 'centername' = OO); geometry:-TangentLine(tgB, B, Elc); geometry:-TangentLine(tgC, C, Elc); dr := geometry:-draw([D(color = blue), lp(color = red), Elc(color = green), A, B, C, T(color = black), H, tgB, tgC], printtext = true); plots:-display([dr], axes = none, scaling = constrained) end proc;
iframes := 10;

display([seq(Fig((1/12)*Pi+i/(10*iframes)), i = 1 .. iframes)], insequence, scaling = constrained);
How to improve this code ? Thank you.

Hi,

Could anybody help me with how to solve the linear time-varying differential equation set in Maple?

A good example will works well.

I've finally decided to see what I could to do with DocumentTools.
I started with these simple lines from ?DocumentTools:-Do 

restart:
with(DocumentTools):
Do(int(%expression, x=%from::integer..%to::integer));
Error, (in DocumentTools:-GetProperty) Attempted to retrieve property of unknown component expression.

These next ones return the expected plot but an alert window 

restart:
with(DocumentTools):
Do(%Plot1=plot(x^2, x=0..1));


However the Click and Drag icon is disabled in the plot toolbar and invisible in the plot menu.

The behaviour of these commands do not depend on whether they are written in a worksheet, a section or a code edit region.

Could you detail me the very first steps do use DocumentTools features?

TIA

DT-Do.mw

Suppose we want to compute the following derivate:

x := [x__1, x__2, x__3]

f := [x__1]

Gradf := Matrix(1, 3, (i,j) -> diff(f[i], x[j]));

Gradf2 := Matrix(3, 3, (i,j) -> diff(Gradf[i], x[j]));

 

my question is how to convert Gradf from [1  0 0] to [1,0,0]?

 

Hello there, 

Would you allow me to ask this question?

What I tried in the worksheet below is to map an expression to x-axis and another expression to y-axis, while both expressions have the same input (or independent variable). However, the resulting plot still shows those two expressions as two separate curves. Would you tell me how to map or assign the expressions to each axis?


 

restart;

with(LinearAlgebra):

Digits := 6:

j := I:

L__adpu := 1.66;
L__aqpu := 1.61;
L__lpu := 0.15;
R__apu := 0.003;
E__tpu := 1.0;
MVA := 555;

1.66

 

1.61

 

.15

 

0.3e-2

 

1.0

 

555

(1)

A__sat := 0.0312500:
B__sat := 6.93147:

I__tpu_compare := (P__t - j* Q__t) / E__tpu;

-(1.00000*I)*Q__t+1.00000*P__t

(2)

I__tpu_a := abs(I__tpu_compare);

abs(-(1.00000*I)*Q__t+1.00000*P__t)

(3)

Phi := arccos(P__t / (I__tpu_a * E__tpu));

arccos(1.00000*P__t/abs(-(1.00000*I)*Q__t+1.00000*P__t))

(4)

E__a__tilde := E__tpu + (R__apu + j * L__lpu) * I__tpu_compare;

1.0+(0.3e-2+.15*I)*(-(1.00000*I)*Q__t+1.00000*P__t)

(5)

psi__at := abs(E__a__tilde);

abs(1.0+(0.3e-2+.15*I)*(-(1.00000*I)*Q__t+1.00000*P__t))

(6)

psi__I := A__sat * exp(B__sat * (psi__at - 0.8));

0.312500e-1*exp(6.93147*abs(1.0+(0.3e-2+.15*I)*(-(1.00000*I)*Q__t+1.00000*P__t))-5.54518)

(7)

K__sd := psi__at / (psi__at + psi__I);

abs(1.0+(0.3e-2+.15*I)*(-(1.00000*I)*Q__t+1.00000*P__t))/(abs(1.0+(0.3e-2+.15*I)*(-(1.00000*I)*Q__t+1.00000*P__t))+0.312500e-1*exp(6.93147*abs(1.0+(0.3e-2+.15*I)*(-(1.00000*I)*Q__t+1.00000*P__t))-5.54518))

(8)

K__sq := K__sd;

abs(1.0+(0.3e-2+.15*I)*(-(1.00000*I)*Q__t+1.00000*P__t))/(abs(1.0+(0.3e-2+.15*I)*(-(1.00000*I)*Q__t+1.00000*P__t))+0.312500e-1*exp(6.93147*abs(1.0+(0.3e-2+.15*I)*(-(1.00000*I)*Q__t+1.00000*P__t))-5.54518))

(9)

X__ad := K__sd * L__adpu:

X__d := X__ad + L__lpu:

X__aq := K__sq * L__aqpu:

X__q := X__aq + L__lpu:

delta__i := arctan((X__q * I__tpu_a * cos(Phi) - R__apu * I__tpu_a * sin(Phi))/(E__tpu + R__apu * I__tpu_a * cos(Phi) + X__q * I__tpu_a * sin(Phi))):

e__q := E__tpu * cos(delta__i):

i__d := I__tpu_a * sin(delta__i + Phi):

i__q := I__tpu_a * cos(delta__i + Phi):

i__fd := (e__q + R__apu * i__q + X__d * i__d)/(X__ad):

i_armature := sqrt(i__d^2 + i__q^2):

i_armature2 := subs(P__t = 0, i_armature):

i_armature3 := subs(P__t = 0.5, i_armature):

i_armature4 := subs(P__t = 1.0, i_armature):

plotA := plot(i_armature2, Q__t = -0.6..0.6):

plotB := plot(i_armature3, Q__t = -0.6..0.6):

plotC := plot(i_armature4, Q__t = -0.6..0.6):

plots:-display([plotA, plotB, plotC]);

 

i__fd2 := subs(P__t = 0, i__fd):

i__fd3 := subs(P__t = 0.5, i__fd):

i__fd4 := subs(P__t = 1.0, i__fd):

 

plotD := plot(i__fd2, Q__t = -0.6..0.6):

plotE := plot(i__fd3, Q__t = -0.6..0.6):

plotF := plot(i__fd4, Q__t = -0.6..0.6):

plots:-display([plotD, plotE, plotF]);

 

plot([i__fd2, i_armature2, Q__t = -0.6..0.6]); #failed, two expressions are plotted seperately.

 

 

 


Best Regards, 

In Kwon Park 

Download Q20210416.mw

What is the reason for the warning?

restart;
f := (x, y, z) -> 1 - x^4 + ln(x^4) + z^4 - y^4 - y^4*(ln(z^4) - ln(y^4));


with(Optimization);
Minimize(f(x, y, z), {x^2 + y^2 - z^2 - 1 = 0, 0 <= x - 1, 0 <= y - x, 0 <= z - y});
Warning, no iterations performed as initial point satisfies first-order conditions
                 [0., [x = 1., y = 1., z = 1.]]

Suppose we have a state vector x=[x_1,y_1,theta_1,x_2,y_2,theta_2] and a vector field f=[cos(theta_1),sin(theta_1),0,0,0,0]

How we can define a vector field with 6 dimentions?

with(VectorCalculus);
f := VectorField([cos(theta__1),sin(theta__1),0,0,0,0],'cartesian ??? '[x__1,y__1,z__1,x__2,y__2,z__2]);

I need to compute the Jacobian of f wrt x:

Jacobian(f, [x_1,y_1,theta_1,x_2,y_2,theta_2]

Dear all

How can I plot vector using maple

Maple return the following error

Error, invalid neutral operator
plotvector.mw

Thank you

In trying to create a plot that shades the region between two concentric circles, I decided to try using the Picture procedure found here: https://www.mapleprimes.com/posts/145922-Perimeter-Area-And-Visualization-Of-A-Plane-Figure-

This is as far as I got.


 

"Picture := proc (L, C, N::posint := 100, Boundary::list := [linestyle = 1])    local i, var, var1, var2, e, e1, e2, P, Q, h;    global Border;    for i to nops(L) do    if type(L[i], listlist(algebraic)) then P[i] := op(L[i]) else    var := lhs(L[i, 2]); var1 := lhs(rhs(L[i, 2])); var2 := rhs(rhs(L[i, 2])); h := (var2-var1)/N;    if type(L[i, 1], algebraic) then e := L[i, 1];    if nops(L[i]) = 3 then P[i] := seq(subs(var = var1+h*i, [e*cos(var), e*sin(var)]), i = 0 .. N) else    P[i] := seq([var1+h*i, subs(var = var1+h*i, e)], i = 0 .. N) end if else    e1 := L[i, 1, 1]; e2 := L[i, 1, 2]; P[i] := seq(subs(var = var1+h*i, [e1, e2]), i = 0 .. N) end if end if    end do;    Q := [seq(P[i], i = 1 .. nops(L))];    Border := plottools[curve]([op(Q), Q[1]], op(Boundary));    [plottools[polygon](Q, C), Border];    end proc:"

L := [[3, t = 0 .. 2*Pi, polar], [2, t = 0 .. 2*Pi, polar]]; plots[display](Picture(L, color = yellow, [color = brown, thickness = 3]), scaling = constrained, view = [-4 .. 4, -4 .. 4])

 

NULL


Clearly, this is not what I want. I am probably missing something simple here. How can I get only the region between the two circles to be shaded?

Download shading_annular_regions.mw

 

 

 

First 434 435 436 437 438 439 440 Last Page 436 of 2241