On the solar differential rotation

(1)

ON THE SOLAR DIFFERENTIAL ROTATION

JUU

\96i

by

Takeo Sakurai

(2)

ON THE SOLAR DIFFERENTIAL ROTATION

by

Takeo Sakurai*

':c Professor, Department of Aeronautical Engineering,

Kyoto University, Japan, temporarilyon leave at UTIAS

Manuscript received November 1965

(3)

ACKNOWLEDGEMENT

The author is very grateful to the National Research Council of Canada for providing him the opportunity to stay at the Institute for Aerospace Studies of the University of Toronto as a grantee of a Postdoctorate Fellowship on leave from Kyoto University, Japan.

Thanks are due to Dr. G. N. Patterson for making the author' s stay at the Institute both pleasant and fruitful.

Finally, great thanks are due to Dr. J. H. de Leeuw for his stimulating discussions and for his friendly encouragement throughout this study.

This paper is dedicated to the late Professor Susumu Tomotika

(4)

'

.

SUMMARY

The purpose of this paper is to re- examine the assumption of

non-isotropic eddy viscosity (by Biermann and Kippenhahn) and to provide a

direct quantitative comparison of theory and observation in thé problem of the

solar differential rotation. The basic equations are taken from the results of

the mixing. length theory of the turbulence and the boundary conditions are

al-most the same as those of Kippenhahn. The important differences are as

follows: (1) We·take {nto acèount the effect of the variation of the density on

the eddy viscosity by assuming that the turbulent velocity is constant and that

the ~ixin-g'length is equal to the scale height. (2) We assumed that the

ro-tational speed on the inner boundary of the hydrogen convection layer is that

of the solid body rotation. This assumption is necessary for our approximate

method since our mathematical problem will become under determinate with-out it. This will also purify our problem by ascribing the reason of the

dif-ferential rotation wholly to the effect of the non-isotropy of the eddy viscosHy.

(3) We applied a perturbation method with respect to the small parameter E

which is the ratio of the thic.kness of the hydrogen convection layer to the solar

radius. This method gives us room to take ipto account the fact that the

ro-tational velocity of the sun is very fast with respect to the speed of the

meridi-onal current. The solution for the first order approximation is obtained by

the usual method of series expansion and the numerical calculations are per-formed on an IBM 7094 electronic computer. The comparison of theory and

observation shows that the agreement is almost complete up to about 40 degrees latitude and the largest discrepancy is 13.4 per cent at the pole,

provided that the thickness of the hydrogen convection layer is 12. 5 per cent

of the solar radius and that the parameter corresponding to the non- isotropy

of the eddy viscosity is 1. 5. Another interesting feature of our results is that

the contours of constant angular veloc.ity in the meridional plane have a buIging

part at high latitude. The distribution of the meridional current is also

(5)

TABLE OF CONTENTS

Page SYMBOLS

1. INTRODUCTION 1

'2. BASIC EQUATIONS 3

3. APPROXIMA TE SOL UTI0N OF THE BASIC 11

EQUATIONS

.

4. NUMERICAL RESULTS AND DISCUSSIONS 21

REFERENCES 26 APPENDIX 1 27 APPENDIX 2 28 TABLES FIGURES 1 to 10 'l iv

(6)

'.

N ote adde d in proof:

On the basis of quite a good agreement of the theoretical and the observational results, the author believes that the solar differential rotation can be explained by the slight non-isotropy of the eddy viscosity. From the viewpoint of the theoretical consistency, however, our theory has a fatal defect: We assumed that the shear stresses vanish at the interface (Eq. (18) ), and also that the rotational angular velocity is constant at the same interface (Eq. (20) ).

The first condition is based on the fact that the molecular viscosity in the

thermally stabIe inner region can be neglected with respect to the eddy viscosity in the outer hydrogen convection layer. In another word, the inner stabIe region

is treated as a sphere of the inviscid fluid. In effect, there is the discontinuity in the meridional current at the interface, since the boundary value of

Us

does not vanish in our result in contrast with the very small von Zeipel's flow in

the stabIe region. There may be the similar discontinuity in the rotational angular velo city. Then, what is the physical foundation of the seéond boundary condition above? This is the question which the author worried about all along his way home to Japan (he left Toronto on the 20th of September in 1965). After he arrived in Japan, he was involved in the problem to get rid of the above

second condition and succeeded in early December in 1965. He gave a lecture on the 21st of December in 1965 at Research Institute for Mathematical Sciences of Kyoto University about his work on the solar differential rotation inc1uding the revised version which is free from the above second condition. He is now writing a paper about his new treatment which will be published elsewhere. In this short note, we want to give a quick survey of our new treatment.

Now, let us consider the second order approximation of Eq. (28)

with respect to E. :

[

~ ~

d

U. -,

1 [ "-

U

;7ï

"leg

~~.

-t-

rr

'Q::

-t ( ; -

~s)

U"j-l-(ZJI;;:.+j,)

l~ ~.

-I-

(J -

25)

iJ~

.

+

s

i&

~

[5iI'/3&

2- {

U&]?7

Si"a8

'd

e

d

ff

Shl

P

J

[

f

~{ (~

a.

U()";,

+

Cosf)

lh \

D: ( ;)

~f)

co.s {)

Lr )

'd

5

iA ()

'(il)

+- "

d

t)

T

~}A

B

'fIJo/

(N1)

"""

Integrating this equation from 0 to 1 with respect to /7. and taking into account the boundary condition (19)~ we obtain the following equation:

(7)

f{'

~-""N

('&0;0

cosB

1Y);p~ L~a"

co.sea·)

:::. /7,

Po

VBu

?;)

e

1- f 111

~

VI'f)j'/o

UlO

C

îJ

(f

+s;;;P

'Ij

\1

(N2)

+({

a~ (~~(+

i1JJdll

By virtue of Eqs. (44) and (45) and the boundary conditions (48)-(51), the in-tegrand of the right hand side of the above equation is uniquely determined as the functional of the zeroth order approximation ~p, Eq. (N2) is, therefore,

taken as the condition to determine the zeroth order approximation without any

ambiguity. Moreover , as is shown by its course of derivation, Eq. (N2) ascertains the mutual consistency of the basic equations and the boundary con-ditions. So that, we call Eq. (N2) the compatibility condition in the second

order approximation. By the similar procedure, we can obtain the compatibility conditions of the higher order approximation. These compatibility conditions take place of the above second condition (i. e. the boundary condition that the

rotational angular velocity is constant at the interface) and make our

mathemat-ical formulation self consistent.

The solution was obtained by expanding each order approximation with respect to the degree of the non-isotropy of the eddy viscosity:

ó=

$-1

.

The solution was obtained up to the order of

d

and

f=

for

0;,

and up to the

order of

b

for ~. The numerical results are obtained for the case with

ri = 0.8 Rj!) and

5

= 0.8, 1. 0 and 1. 2 on FACOM electronic computer. This

value of r i is chosen on account of the theoretical results on the hydrogen con-vection layer by Schwarzschild, Howard and Harm. ':< From the results shown in Fig. 1, we can conclude that the solar differential rotation can be explained by the slight non-isotropy of the eddy viscosity in the hydrogen convection layer.

We can also conclude that our results based on the application of the above second condition can safely be applied to the astrophysical problem as far as we take into account the inconsistency in the formulation.

i,< M. Schwarzschild, R. Howard and R. Harm : Astrophys. J., Vol. 125,

1957, p.233.

(8)

Rotational Angular Velocity

1·5

1·0

0.5 C

.-v

r;=

1.0 L---:::::

,... ~ '(,: O.

8=0.8

8=

'.0

~

-

_Ir:

_Q_.

-

_-

I- _

_r.

₌₌_/.0

--

Ob

s;;'vatIön

-g=/,';.

30

0

60

0 _~o Latitude

FIG. 1 NUMERICAL RESULTS OF THE ROTATIONAL ANGULAR

VELOCITY VS LATITUDE AND THEIR COMPARISON WITH

(9)

=

f

S

®

't

ri-re

R· ₁ = € == ~o ri. ho

=

SYMBOLS

system of spherical polar coordinates velocity

density pressure

gravitational potential stream function

viscous force due to the Reynolds stress, Eqs. (7) (9)

Reynolds stress, Eqs. (la, 1) - (14)

eddy viscosity corresponding to the vertical

component of the turbulent velocity, Eq. (15)

the ratio of the horizontal part of eddy viscosity

to the vertical part, Eq. (16)

turbulent velocity the mixing length

components of the diagonal tensor corresponding

to the turbulent pressure, c.

f.

the expressions

in the text.

ratio of specific heats

inner radius of the hydrogen convection layer radius of the sun

n

/re

small parameter used as an expansion parameter

in our approximation method, Eq. (23)

ri - re

=

thickness of the hydrogen convection

layer, Eq. (23)

(10)

/

suffixes 0

tildes

G

₌

M

₌

used for the reference value of the quantities

used to indicate that the quantities are non-dimensional Eq (29)

J2-

=

Po

ho

)0 U_r;

ri.

Eq. (30)

J

3

)0

Po

=

po

J

q

Áo

_L{rlfo = _Eq. (31)

r"

iAto

Rr

₌ -Á.2. C>

I~

r

"'Zo

fAro ...,

A:

Re

= Eq. (32) 2. ~o

re

~ro

R~

-

it~

R

1> 10

rL

U</J_o

1

0 =

f"

J,."

t..f

"'0*

Eq. (33) gravitational constant total mas s of the sun

..

(11)

suffixes 0, 1, 2 et ceterj

r

used to indicate the coefficients of expansion of the non-dimensional = lformula

quantities

[

('(-/) G M

fo

Jto}1=ï

j

= 5"

Y

po

rt

C

I +

E:)

"...,

(0-1)

J=r

J

_5"

=

1+ E.. ,...., _,...

r,

=

1-

r

"-"

F / (~)

_-

integration constant of Eq.

~

-

-'--

L

;tt)

=

r,

'a'-I d ( ... '1"-1

- r ,

_dr,

L

J ::

2 J

5"

J

,2. t: 2-6

:95"

.P'f

J

7 = 2

J

S"

J"f

J5"

j

=

1r(1-2S) 8

_(t-/)

J

₃

(

)2.

2.

2-.

-

'+

I

+

E.

U.

~o

({

"-

/)2. /..,( 2.

r*o A

.

-

J

₇ cos 9 B

₌

Jv

sin 9

c

:

-

./)6

cos 9 x (67) Eq.(.56) Eq. (64) Eq. (60) Eq. (69) Eq. (70) Eq. (71) Eq. (72) Eq. (77), Eq. (84) Eq. (85) ~q. (86)

(12)

D C sin 9 Eq .. (87) ,..., E

-

C

F,

(9) Eq. (88)

Ó

n. m

-

Kronecker' s delta ~

TL -

e

Eq.(111)

G=

::: latitude measured from the equator

(13)

According to the estimates of Cowling (Ref. 4) the lifetime of the inequalities of the rotational velocity in the thermally stable region of the sun is of the order of 10 11 years, provided that the temperature is 10 5 degrees and the density is 10- 4 gm/ cm 3 . Since the age of the solar system is considered to be of the order of 10 9 years, the rotational motion of the stable region must

be a relic of the past which indicates its state when the solar system was first

formed, as Cowling suggested.

In the thermally unstable hydrogen convection.layer outside the

stable region, on the other hand, the corresponding lifetime based on a mean

turbulent velocity of 10 5 cm/ sec and a mixing length of 10 9 cm is of the order

of 10 7 years. Hence,the inequalities of the rotational velocity w.ill be redis-tributed by the effect of the eddy viscosity in the course of time since the solar

system was first formed. Moreover, H. W. Babcock (Ref. 2) estimated this

lifetime as equal to a few thousand years on the basis of his theory of sunspots,

the estimate of the energy Of the solar differential rotation (10 39 ergs) and of

the magnetic energy required in one sunspot cyc1e (10 36 ergs). There is

another estimate by Ward which shows this lifetime to be only about 100 days (Ref. 8).

We can therefore safely conclude that the rotational motion of the outer hydrogen convection layer must already be at its final asymptotic

state while that of the inner stable region is still in its initial statë.

Conventional hydrodynamics based on an isotropic viscosity shows that the asymptotic state of the rotational motion must be a solid body rotation. The rotational motion of the hydrogen convection layer, however,

shows a considerable deviation from the state of solidy body rotation as is

shown by Newton' s formula of the solar differential rotation:

(1)

where uJ is the angular velocity of the rotational motion of the sun and

p

the latitude.

To explain this characteristic of the rotational motion of the sun, one approach is to assume non-isotropy for the viscosity. This may be

based on. the following assumptions:

(1) The fluid motion in the hydrogen convection layer is completely turbulent

(2) The turbulent eddy viscosity, which is expressed by the formulae of the

mixing length theory of the turbulence, is many orders of magnitude larger than the usual molecular viscosity.

(14)

(3) There must be an inherent difÎerence between the ve.rtical and the

hori-zontal turbulent motion because of the 11.0n- isotropY in the body forces

such as the gravitational force and the force by the magnetic field.

Such an approach was initiated by Biermann (Ref. 3) who showed

that although the steady state solution of the rotational motion of a fluid sphere

with a given angular momentum is a solid body rotation in the case of isotropic

viscosity, the solution for the case of non- isotropic viscosity is much more

complicated. Kippenhahn further developed this approach by investigating the

case with more realistic boundary conditions and concluded that the theory of

non- isotropic viscosity gives solutions with features characteristic of the solar

differential rotation. While these theories are very instructive, both have failed to give direct quantitative comparison between the theoretical results

and observations. Since the assumption of non- isotropic viscosity seems. very

promising to the autho~, it is very important to determine its applicability

from the quantitative point of view.

In this paper, we win re- examine the problem, which was

proposed by Kippenhahn by using a little more realistic assumption on the eddy viscosity and a completely different method of solution.

As in the paper by Kippenhahn, we restrict ourselves to the

consideration of the axisymmetric fluid motion in the th in hydrogen convection

layer. The boundaries

ot

the'layer are considered to be two spheres with

radii ri., and re. (re < re. ), _{respectively. We take Kippenhahn ' s boundary}

conditions of vanishing shearing stresses. We also take a phenomenological

point of view without any regard to the detailed mechanism of how the non~

isotropy of the viscosity is sustained. Another assumption which we borrowed from Kippenhahn is that the convective layer is adiabatic and that we can

neglect the components of the Reynolds stress which do not depend on t!1e mixing length. However, we purify the problem by assuming a solid body rotation in the inner region. This condition is necessary for our method of

approximation since our mathematical problem win become under determinate

without this condition. This corresponds to attributing the reason for the

differential rotation wholly to the redistributing effect by non- isotropic viscosity

rather than to the complicated unknown processes when the solar system was ""

first formed. Secondly, we assume that the mixing length of the turbulence

is equal to the scale height, which seems a reasonably weU estahlished assumption from the astronomical point of view. We also assume that the turbulent velocity and the ratio of the vertical and the horizontal component of the turbulent velocity are constant. FinaUy, we use the first approximation of a method of expansion in terms of the ratio of the thickness of the hydrogen

convection layer to the radius of the sun. 'This final point makes the analysis

completely different from that by Kippenhahn and win clarify the reason why

(15)

2. BASIC -EQUATIONS

Let us consider the turbulent motion of a compressible fluid in

a spherical shell and assume that the averaged quantities, such as averaged velocity and density, are axisymmetric and steady. The basic equations for

the averaged quantities are obtained by the mixing length theory and can be

expressed,in terms of spherical polar coordinates, as follows:

cd

t..i_r _::

U!

_f

L.i;

-'

dP

_.&

+

ik

tAt

_r

r

;:)r

'dr

f

(2)

d

l.i~ _L.j~

_cos

_e

_(,.,(r _~&

_L-

_d_P

Ba-- = -f

clt

_r

sIn

e

r

fr

~9

fr

(3) (4)

Yr=

I

ó-!P

"'1& ::.

,

èJ!p

fr2.

sIn

e

de ) _frsl~e _~ ( 5) where

ft=

LA,.

..L.

₊

~

L-dt"

_r

_'de

( 6)

and (r)&.(j'J) ) (L4,.Il..l&,

tA.)ol ),

p)

p,

Qn~

cl

aretheco

-ordinates, the velocity, the density, the pressure, the gravitational potential,

and the stream function, respectively. The quantities

R,.., Re

and R~ represent the Reynolds stress and can be expressed as follows:

(7)

o

(8)

(9)

(16)

where ...,1"

- 21

d /..ft"

Cr

=

(10. 1) ~r

...

~

1:/'

- 1l.

[:r

(r~s)

+

S

~~c

- 2.

S

(..j. ]

Cr

=. ~= (10) r~

r.r

Trp

::'

t.

=

-

ra..

s :'" t

e [;

r (

r

S J r, 6

~,..)

- 2 S

SIn

IJ

~

} ( 11)

r

_r'

SI""

e

(12) (13) (14)

are the components of the Reynolds stress and

(15)

s

=

(16)

are the eddy viscosity corresponding to the vertical component of the turbulent

velocity and the ratio of the horizontal part of eddy viscosity to the vertical

part, respectively. Here, (IA,..,. J (..( • • ~

1..1,. )

and

eB)

are the components of the

turbulent velücity and the mixing length, respectively, and the assumption is made that there is no preferential direction in the horizontal components of

turbulent mot ion. The equality of the vertical and the horizontal components

of the turbulent velocity leads to the case with S

=

I ,

that is to say the case

with isotropic eddy viscosity. Hence, the assumption of non-isotropic eddy viscosity relies implicitly on the assumption of the non-isotropy of the turbulent motion which may be expected from the non -isotropy of the body

(17)

....

.'

forces such as the gravitational force. These expressions are borrowed from those in the book by Wasiutynski (Ref. 9) and, as in the case of Kippenhahn,

we have neglected the effect of ordinary molecular viscosity and those terms which do not depend on the mixing lepgth: i. e. i'rrr

=

f

'-1;.

'

11.₈ "

r-a.

f<J; ...

and ~~

=

ra. sIn"

e

fC,..fîlt

Sinçe the mean thermal energy per unit mas-s

is estimated tobe of the order of 10 13 to 10 14 erg/gm while the mean kinetic energy of turbulence is about 1010 er:g / gm in the hydrogen convection layer,

the omission of terms like 1rrr is well established. Bècause otH-1ilits,,·simplifio.$.Hon

we can discuss the average velocity field separately without any regard to the kinetic energy of the turbulence itself. Of course, for the analysis of the velocity field~the reasonable assumption is made that the turbulent shear can

be expressed in terms of a mixing length and that the turbulent velocity is constant. This situation is very fortunate for our approach since, in the ex -perience of hydrodynamicists, the mixing length theory is a very good ap-proximation when we apply it to the problem of investigating the average veloc ity field only.

In place of the energy equation, we assume that the fluid is adiabatic in the hydrogen convection layer. :

(17)

where

'ó

is the ratio of specific heats.

Since the eddy viscosity in the outer layer is many orders of magnitude larger than the molecular viscosity prevailing in the inner stabIe region, we can approximate the condition of continuity of the shearing stresses by the vanishing of these of the outer layer on the interface with the inner

sphere:

1:;

= 0 on

r

=- r~ (18)

The same conditions must also be applied on the outer edge of the layer since in our analysis, the existence of the tenuous chromospherical

and coronal medium is neglected. The meridional current in the stabIe region is estimated to be quite slow (of the order of 10- 8 cm/ sec), so that we can assume that the meridional current in the outer layer cannot cross the interface .{as well as the outer boundary~. :

p=

0

ol') ₍₁₉₎

These two boundary conditions and the assumption of the adiabatic relation are borrowed from the paper by Kippenhahn since they seem to be quite reasonable.

(18)

These two conditions are, however, not enough to determine the flow field of the hydrogen convection layer at least in our approximation method. Hence, we assume that the distribution of the rotational velocity on the inner boundary is equal to that of a solid body rotation. :

(?O

)

This assumption stems also fr om the idea to purify the problem by ascribing the reason of the solar differential rotation wholly to the effect of the non-isotropy of the eddy viscosity rather than to the complicated un-known processes which occurred when the solar system was first formed. This corresponds also to an interpretion of the pr~sent state of rotational motion in terms of the final state in which the inner stabIe region itself rotates as the solid body, provided that the meridional current in the inner region may still be neglected.

height:

Finally, we assume that the mixing length is equal to the scale

®

oe

f

dP/

/d,..

(21)

and that the turbulent velocity is constant over the entire hydrogen convection

layer:

=

constant (22)

The first of the above assumptions is an extentipn of the generally accepted assumption on the size of granulation to our problem.

Applying the ordinary transformation used to non-dimension-alize the basic equation:

r = r{. (ITé

r )

[

/'0

t

=

re

- r;:

(23)

r:t,

' V ure> ,-v N

u

r Z

u

to ur )

u

_e tie

Urp

=

u,po Uf/J (24)

E.

(25)

where suffixes 0 indicate the reference values of the quantities and tildes

(19)

the following equations:

_ -

J~

&

[0

t f:

~)'

(t

Ü ;.

J,

~

f )

1

é.'i.

~,

.L f

sinect

;!

+

;~

d

W)}

sine dG

L

~e ~

+

J~ ~

r {(

I

+

é

r)i.

~}

+

I

Ü

~ L(,. -..

}

I+Er 9 ~

d

f

J'

coS

e

U,.,2. -

E-

U,.

iJ

e }

o-r

l

I SIne .,..

f

sIn

e

'7 (26)

,

(27) (28)

(20)

where use is made of the divergence and the

<P -

component of the rotation of the equations of motion instead of Eq. (2) and (3) and

J,

=

.Ir..

t,A<Po (29)

rt. U

ro

J

=

pfJ

. .J....

J

3

fo §o

20

=

₍₃₀₎_.

2-fo

Url)

r·

I.

_Po

J

:=

Jr.~

fA,.,.o

(31) ~ 1...1

"0

rt.

J..:

RI:

tf'-

lt:

Re

,....-

~: R~

- . "

R4>=-Rr=-

Re

:: ₍₃₂₎ "'(. fA ....

1-

"t

_(.;(ro ")

"1-

ri.

{.)4>0 (33) Here use is also made of the assumption in. Eqs. (21) and (22).

Now, let us suppose a situation in which the right hand side of Eq. (28) is predominant, say for example of the order of 10. Then, (.;("'0 must be of the order of 103 cm/ sec, since /...Ir.o is of the order of 104 cm/ sec from the results of observation and

.lt.

is of the order of 10- 1 of the solar radius.

This leads to the result that the predominant term in Eq. (27), on the other

hand, is the third line of the left hand side, since L,(~o is of the order of 10 5 cm/ sec and hence

J,L

is of the order of 10 2 . Hence, even if we can assume the predominance of the viscous term in Eq. (28), the same is not at all true in Eq. (27). This situation does not essentially change within the reasonable assumption on

t..1".0 ) ..1,,0 )

e. t. c.. It is the reason why

Kippenhahn' s method of approximation failed for the case of the sun, since his method relies on the implicit assumption of the predominance of the viscous

term also in Eq. (27). The same estimate leads to the overwhelming

pre-dominance of the first line on the right hand side of Eq. (26) on the basis of the fact that

P:!fo

is of the order of 10 13 cm/ sec and hence that

J2.,

is of the order of 10 9 .

From the above estimate, we can safely neglect the effect of

fluid motion for th~ pressure. Hence, the basic equations to determine the

pressure distribution is completely equivalent to those in the theory of the

structure of the sun in the case of convective equilibrium. Since the mass in

the hydrogen convection layer is many orders of magnitude smaller than the

total mass of the sun, these equations are, as is shown by Str~mgren (Ref. 6)

(21)

(34)

(35)

where G and M a r e the gravitational constant and the total mass of the

sun, respectively, and use is again made of the original, dimensional ex-pressions for the physical quantities.

Before proceeding further into the approximation of Eqs. (28) and (27), note that the thickness of the hydrogen convection layer is about

10 to 20 per cent of the solar radius and that

e...

can hence be considered as

a small parameter. On the basis of the rough estimate above, let us assume the orders of magnitude of parameters as follows:

By applying the usual procedure of a perturbation method on the non-dimensional quantities: -v

---

"'-U

,..

= vr

ro

+

E.

Ur

,

+

-.

.

~

.

" '

-a •.

....

~è

=

+

€.

Lle, ....

.

,.

..

.

-

_crct>c

₊

...,

U

fI = E

LIep,

+ ....

#

-

...,

-r

=

fo

+ é.

)1

+ .

t • , ,.

--

~

f,

~

...

~

=

~o

+

..

9 (36) (37) (38) (39) (40) (41)

(22)

we can obtain the basi~ equations of the zeroth and the first order of approxi

-mation:

zeroth order approximation:

(42)

(43) fir st order approximation:

.

l{~

d

ciA ••

U;.) _

'iJ

d

lA,.}

+-

~;)

_

{k

~

(:n.

d

Kro'}

=-

0

2. I SIr)

e

~ ~

••

d

e

~ d"

r.

d

r

C' ~

r ')

(44)

-

-' V

d

U~.

+

"" d

t.i

te.

-

" - _cOS'_S

U

rQ

U.

co

+

U

e• ~,. S"'ÏrÏë

~r. 6~

J"I-

~

[;.

[

~~~

+

(I -

2 5)

i:T,.} ]

-

_r

_àr

(45)

The same procedure leads to the approximate version of the boundary con-ditions as follows: (46) (47) (48) (49) (50) (51)

.

(23)

-As the boundary condition for Eqs. (34) and (35), we assumethat the pressure

and hence the density vanish at the outer-edge of the hydrogen convection

layer:

"

-p=o

°11 Y'=

(52) This corresponds to neglecting the existence of the tenuous chromospheric and coronal medium. From this condition, we need not impose the boundary conditions related to the vanish of the shearing stresses on the outer edge. The only necessary boundary conditions are as follows:

(T;

co

~o

,....

on r::: I

3. APPROXIMATE SOLUTION OF THE BASIC EQUATIONS

The solution of Eqs. (34) and (35) satisfying the boundary condition (52) is easily obtained to give the following result:

(54)

"-:: J

!, -

~

)t=i

f

'" ,-,

+

€.

r

(

5 5 )

which is also written in the following form:

(57) ...L-... A - '(-,

fel

=

o::P~

r,

(58) (59) where (60) 11

(24)

By substituting these into Eqs. (15), (21), (22) and (41),we obtain the foilowing expression:

-

X--

Jr

...

'6-1

1

0

~

-

'(- ~

f. (

I-~)

"7,

₌

_{'Ir - ,}

J~

=

(i-/)J.

_,+~ (61) (62) (63) (64)

The solution of the basic equations of the zeroth order of approximation Eqs. (42) and (43) satisfying the boundary conditions (46) and (47) is easily

obtain~d to give the following result:

(65)

The angular velocity distribution. in the zeroth order of

ap

-proximati0n 0n each sphere corresponding to each value of ,.., , is therefore

like that 0f the solid body rotation as in the case of Kippenhan who· took the

f01l0Wing expression as his zeroth order of approximation:

~

)-2.<

1-5)

lA.o.:" (,

+

e.

F

sin

9 (66)

The important difference on, the other hand , Ïa:thm GUI<~efl.Oth

order approximation is constant with respect to

r

while his inc1udes terms

depending on

F .

This sterns from the fact that the centrifugal force term

is predominant in' Eq. (27) rather than the viscous force term.

By substituting Eq. (65) into Eqs. (44) and (45), we obtain the following equations for the first order approximation:

(67)

(25)

(69)

(70)

(71)

(72)

Before going further into the solution of these equations satisfying the boundary conditions (48) - (51), let us consider the relations between these boundary

conditions. By integrating Eq. (68) from 0 to 1 with respect to

r;

af ter the

multiplication by

r

,~I ,we obtain the following equation:

(73)

Using the non-dimensional version of Eq. (5) in the zeroth

order approximation:

,...

p.

= ; ..

;;,,9

r

cr ••

r.t:ï

oI~

(74)

where the integration constant is determined by·the boundary condition (54)

to be zero, and the boundary condition (48), Eq. (73) reduces to the

following:

S In

a

(75)

which is nothing but the boundary condition (50). Hence, the boundary

conditions (48) and (50) are mutually dependent within our first order ap

-proximation. By using the<qt>ë~~~~~~-Eq.(68) and taking into account Eq. (67)

we can obtain the following equation: . .

L L

I:1,

I

= -

J,

c.o

S'è. 9

['LÁ

~,

1-

FE; (e)

+

r,

S I f)

e }

₍₇₆₎

where

J" -

(77)

(26)

This equation and the bounda;ry conditions (48) - (51), (53), and (54) being independent from

J".

and

Js-'

the same is true for the rotational velocity and hence for its angular velocity within the first order approximation. This shows the wide applicability of' our approximation method, since we need not

worry abol.lt the practical choice of the values of these parameters. The

- - A

quantities L.(rrJ and t...1eQ are, on. the other hand, p:r:!?P0rtional to éI~ and in-dependent from j~, since these are related to

tAtP,

througp Eq. (68).

substitution of:

.

into Eqs. (67) and (68) leads to the following:

of which (79A) is equivalent to, the following:

7(6)= ....

.sIne

tca)

= -

F,

Ce)

(78) (79A) . (79B) (79C) (79D)

By (79C),

S

in (79B) must be equal to I so that (78) may be the solution. By the boundary condition (51)~ ~(9) is determined as follows:

fCe)

=

sin

e

Finally, the ·substitution of (79C) and (79E) into (78) gives the solution for the case where the assumption. (78) is admissible, i. e. for the case with

5 ::-

I

.-..;

-L( ro

=-

0 )

_L{&o

=

_{0 ')}

(.../~I -

=

,...,

r-

s, ,..,

I @

The above results are summarized as follows: The meridional current vanishes for the case with isotropic eddy viscosity and, simultaneously, the rotational velocity reduces to that of solid body rotation.

expand î1

'""I f> I

Now, taking into account the botindary condition

_,...

. (53), let Us

.

(27)

(80)

.,., .. 0

(81)

The substitution of these into Eqs. (67) and (68) leads to the following: where A

₌

B

=

c

=

D

₌

E

=

~"") l(~~I)

't -

~}

t - I

&=+

IJ

f~

..

1)0'

-0-z,}

0-1

J

₇co

~

e

Ja

s

i

f1

e

- J

₆

c.~

s

€i

Cs, n

9

-C

F;

(s)

and ~'I'I.)IW\ is Kronecker' s delta.

(82) (83) (84) (85) (86) (87 ) (88)

The solution of Eqs. (82) and (83) is easily obtained to give:

(89)

(90)

(91)

(92)

The substitution of these results into the boundary conditions (48), (49) and (51) gives foUowin,g simultaneous algebraic equations to de -termine the unknown consta.rl.ts

f..,

1

0 and

E :

(28)

(93) (94) "

-A3

fo

+

B3

~()

-ol-

C:; E =-

0

₃ (95) 00 ('6 _ d"L~

11"" C

'YY

A

I ·

I

+

~

l@"",).lr!

(2.

't - ') - - -

{~#I.

Y -- (2"" -1 ) } (96)

..

_ 00 f'{-/)'~-'

e'f'v

C'YI-'

B, -

~, f('2~_t),I},((2r-/)--- [ca~-I)~-('2.~-'2.)J

(97)

Cl

=-

11,-/

c

(98)

P - -

[B

+-

QD

Je

+

12-I -

A

(r -,)

A C i C

(99) _ GOO

Or -

rt

n-, /1"'-'

C,.",

Aa -

ho

{(2,"-C.)n'(('2.~-J)

__ -

fc,,,.,-,)f-

(aM.-'2.)}

(100)

B :

Ë

(o-I/-~ Ii~

c""

2.

~~I

l<2"t.-I)!}"'t(2t-/)---

t?'''Yto-(2'1'1-I)]

(101)

(102)

(

8 )JO)

D-a. ... -

/3a

ft

+ ('(-

f)

AC,

(103)

~

(t

-1

)""fI.-'

11

Hl-I

C

~ -

~

Ü2.If1.-I)!}'7f(21f-I) __ - [2"IV)'-

(2.~-')}

(104)

..,

(r

_,t"-n.

,

p,#f'I.

C ""

B?

::

I

+

~

i""

...,0: I [(2m)

l} '(

ci-r-I) - - - (

(2.11\.+ 1)0 --2~} (105)

C

₃ :=-

A?

C

(106)

D,

-

(B~

-

t)(

t

-+

_('(-/)11

~D

_c

)

( 107)

By substituting the solutions of Eqs. (93) to (95) into (80) and ' V _

(29)

-quantity Urccan also be determined by the substitution of (81) into (74) and the non-dirnensional version of (5). Since the convergence of all of the above

series is very rapid and the expressions for each terms are systematic, the above method is very convenient for computer work. From the point of view of investigating further the details of the flow near the equator, we use also

another method of expansion with respect to the latitude:

(108)

(l09)

(110)

where

è==f-e

( 11)

is the latitude. By substituting these expansions into Eqs. (67) and (68), we obtain.the following equations:

-L

U~

_1'1'\.

=

where A

-r. [

r

2

rI_

J,

I ó"~ - ~ 17 (113) (114)

(30)

The same substitution into the boundary eonditions (48), (49) and (51) leads to the following: (115)

-on

r,

-=

₍₁₁₆₎ (117)

where we assumed Eq. (128) for the derivation of (115).

By equating hl.. equal to zero in the above equations, we ean

obtain the following:

L

LAeoo

= -

J,

(~+ U~'o

.,.

/='0)

(118) (119)

o

(120) "'"

on r,

= (121) (122)

Using the relations derived in the Appendix 1, the solution

of Eq. (119) satisfying the boundary eondition (122) and (53) is obtained as follows:

(123)

Substituting (123) into (118), we obtain the following:

( 124)

Using the same relations in Appendix 1 as used above, E;q. (124) is integrated to give:

(31)

The quantity ft.ois determined by the boundary condition (121) as fo1l0ws:

J,('t-I) _

1 -

[I

+

'oe((-I)

l

:

'(

a'ç-I

0

5

(126)

Substituting (126) into (125) and integrating, the following is obtained:

(127)

Before g~ing furt!!,.er into the determination of ~o , we must

consider the relation of

U

eo to

tflo

.

in our expansion procedure. By the

substitution of the following expansion:

--

-

- ( -

-air:

)

'Po

=

~

P

00

+

e

:t' f) I

+ - - -

::> (128)

(109)and the expansion formula of trigonometrie function into Eq. (74), we

obtain the following relations. :

'.

.

F. -

j

,

~

Ueoo

~rf,

POf)

=

(129) 0

"'-~PQ)Jr.

,..,

r'l,..

-~DI

=

ft)

(U&OI-0 ( 130) ,...

fl

--

Usofr:-j)

Po,."

=

'fo (

D'

tON\. -

+

2.

0

( 131)

The substitution of (127) into (129) leads to the following:

(132)

where use is made of the following relation:

(133)"

The quantity

era

is determined by the boundary condition (120) as follows:

(134)

(32)

Now, Üet>C) and

POtJ

are obtained by the substitution of (134) into (127) and (132):

(135)

(136)

~

lAf'O is also obtained by the substitution 'of (136) into the non-dimensional

ver sion of 5:

(137)

....-where the following assumption on

U"'o

is made in conjunction with (128):

(138)

Taking into account the definition of

Je ,

we can obtain the following

con-c1usion as in the case of Kippenhahn: if the friction is bigger in the radial

direction (0 < S < I) the meridional current rises at the equator while if the

friction is smaller in the radial direction ( S ;> I .> the meridional current faUs

at the equator.

.

By repetition of the same procedure, we can proceed

straight-forwardly into the lligher order approximations. Since' the procedure becomes

much more tedious, however, only the results are tabul~ted in Appendix 2. .

...

"'-By substituting the value of V(I>I/ at

r,

==

0

from the expression in·Appendix 2):

1-

2.

5

2.

•

(which is obtained .

(139)

into Eqs. (108) and (39), the distribution of the rotational angular velocity on the outer edge of the hydrogen convection layer, within this order of ap-proximation, is obtained as foUows:

(33)

This leads to the 'following conc1usion: The rotational angular velocity on the surface of the sun increases with the latitude in the case where the friction. is bigger in the radial direction ( 0 <. S < I), while the situation is reversed in the

case where the friction. is smaller in the radial direction ( I <. $). A similar consideration leads also to the following conc1usion: The rotational angular velocity on the equator of the sun decreases with increasing radius in the case

where the friction is bigger in the radial direction ( 0 <. S

<

I ),

while the reverse is true in the case where the friction is smaller in the radial direction

(

I

<

S ). The first conclusion which seems to be restricted to the

narrow zone near ·the equator because of the situation that we used an approxi-mation valid only for the very low latitudes proves to be true also for high latitudes as will be shown in the next section.

4. NUMERICAL RESULTS AND DISCUSSIONS

Before jumping directly into the performance of the actual

calculations, let us pause here to investigate some properties of the solutions. Equations (94) and (95) can be rearranged as follows:

02.

(141)

(142)

The solutions of these equations are obtained as follows:

(143)

1!J..

B

'3

ti

'ft (3a

C

c,...

êj. -

[~

+

~~/fAC

J[

t

IJ,

B~~

1'1.8,- -

I }

(144) The substitution of (143) and (144) into (93) leads to the following:

(145)

Here

~

,

*" '

~'ë.

and

(B?J-

tl/AC

are the poweJ;. series with respect to

Ifc

starting from zeroth order terms as

C~-~({

,

~-

'Yf '(

(2.r -I)

1 .

('(-'V[rcz.

_{r -}J)] and

tt-

J)/2.[d'(aY-I)(30-2-'1} , respectively.

Taking account of these facts and the following relations:

(34)

AC

:. _ q (

I

+

é ~2. E. Z 2:

U

2-<b0

cos'a e

(J-/) (,{ r~o (146)

.li.+

'(0

_-

E'tUr~o

$

in

e

_(1-

_sj

Ft &-/)/1-C

_U

_+~)

_L.1

_1"0

cos e

(147)

A

2 C

1+

EJ Ur::g

cos

e

Cp

-I)

ê. CA r ~o (148)

we can reach the following conclusions:

...

-(1) The meridional current

U

ro and UaC) is completely determined by the

zeroth order rotational velocity (65) and the boundary conditions (48) and (49)

without any regard to the imposition of the boundary condition (51). Hence.

a slight deviation from the state of rigid body rotation in the inner stabie

region does not change the meridional current within our order of approximation .

. J .. These are proportional to ( 1- S ) and inversely proportional

to Reynolds number ( ; ; ' ) . The distributions of these depend on ifC

which in turn is proportional to the square of the proportional factor of

e..

as

is shown in (146). The qU8ntity

aso

is proportiémal to sin

e

cos

e

which

' V

shows, as is expected, that U&o vanishes on the equator and at the pole.

(2) The rotational velocity U'IH is decomposed into two parts. One is

completely independent from any parameters corresponding to the first term in (145). The other is proportional to ( / - S ) while independent from

Reynolds number. On the outer edge of the sun, the deviation from the state of rigid body rotation is equal to this second part. The velocity distribution

correspondingto the second part depends on

Ifc,

as in the case of the

meri-dional current.

By the above conclusions, we need not scan the velocity

dis-tributions with respect to all the values of S and

J.,.

while it is necessary

to investigate the effect of

E.

on the velocity distributions. We assume as

typical examples, the foUowing values in the actual calculations

5 =

0.5""

and

/.5" ₍₁₄₉₎

(35)

r·

L

o-e, 0-825", 0-85"0) 0-87S", o-!JO) 0925

ot:

"è

(151)

The actual calculations are performed on an IBM 7094 electronic computer at the Institute of Computer Science of University of Toronto and the results are tabulated in the Tables.

The most important feature of observations is the deviation of the solar rotational motion from the state of a solid body rotation. Hence,

from the strictest point of view, it must be this deviation rather than the total

value of the angular velocity that is to be compared with the ob servation. This

point isalso appropriate from the above consideration of the properties of the

solution. Hence, we take the deviation of the angular velocity from the

reference value on the equator both for the theoretical and the observational

results and take the ratio of these. Thi!3 is shown in Fig. 1 where the ratios

of the theoretical values and the observed values are shown. Taking into

account the above conc1usion that the deviation is proportional to s-I , this

ratio must be constant if the agreement is complete. Actual values for these

ratios are positive for the case with

s

> I and negative for the case with 5 , I

showing that the former cases are more realistic. In all cases, the absolute

values decrease rapidly with increasing latitude. This corresponds to the

situation in which the deviation of our solution from the solid body rotation dies away more rapidly in comparison with the observed results. The con-stancy of the ratio becomes better, and simultaneously, the absolute value of the ratio decreases as the thickness of the hydrogen convection.layer becomes thinner.

These features may correspond to the situation in which the

rotational motion approaches to the state of solid body rotation and hence that

the ratio approaches to the constant value of zero as the thickness decreases.

Therefore, the preferred case must be chosen by the criterion that the ratio

is of the order of one as well as relatively constant. From such a· point of

view, the case with RL=' o·E'7S" is the best. This conclusion is plausible since

it is expected by the theory of the stellar structure that the hydrogen con-vection layer of the sun is of ten to twenty percent thickness of the radius of the sun (ReL 6). Another comparison of the theoretical and the observed values are shown in Figs. 2 and 3. In Fig. 2, the ratio of the total angular velocity are shown. Corresponding to the conclusion above, this ratio is

almost one in the case with ~ ~O·B7:; up to about 40 degrees latitude and then

rises gradually to arrive at the value of 1. 134 at 90 degrees latitude. This

figure also shows that the ratios for the case with S;:::./ and s <../ approac~

from opposite sides to the value of solid body rotation (marked S. B. in the figure) corresponding to the vanishing thickness of the layer. Figure 3 shows

the theoretical results for the angular velocities normalized with respect to

the reference value at the equator with the observed values shown by the dotted

line. From this figure, we can conclude that the agreement of the theoretical 23

(36)

result with the observed result with respect to the angular velocity is complete up to about 40 degrees latitude and reasonably good for the higher latitude for the case with S

=

/1 S"" and

Ri.

= ·0,

e

7 S"" •

As was discussed in the irlitial part of this section, the flow

properties change maÏ.i11y with the sign of the parameter S - , , so that the

difference of the flow properties may be illustrated by comparing two cases with S

>

land S <. I. This is done in Fig. 4- 9. Figil.re 4 show s thes"

variation of the angular velocity with latitude with

?; •

as the paramete:r-r~r the case wUh S = I,

'5"

and

RL

= 0·87

s:

The angular velocity decreases with increasing latitude for all

values of

r;

.

correspondjng to equatorial acceleration. while the situation

for the case with S:-

o,!)

is reversed as is shown in Fig. 5. Figure 6 shows

the variation of the angular velo':!ity with

r.

with the latitude as the

para-meter. As was mentioned in thc conclusions related to the expansion method with respect to the latitude. the angular velocity is the smallest on the inner

boundary and increases with radius in the equatorial region. As the latitude

increases, however, this situ:?tiol1 gradually changes and at the latitude of

ab'~ut 40 degrees, the gradient of the angular velocity with respect to the radius becomes zero on the surface of the sun. Proceeding further into the

higher latitudes, there appears a maximum in the curves of the angular velocity. The meaning of this feature wiU become clear by the contours of

the constant angular velocities in Fig. 8. Up to about 40 degrees latitude. the contours have relatively simple monoton:i.c shape as was assumed by

Babcock, .(Ref. 2). Higher than this latitude. however, t11e contour s have

buIging parts. The contour corresponding to the horizontal tangent of the curve (with 90 degrees latitude as the parameter) defines the boundary of the

high latitude region where the angular velocity decreases with increasing

radius. If a sun spot can be considered él,S the portion of the horizontal mag-netic lines of force which is submerged at a shallow depth Ln the hydrogen

convection layer and thrown up by the effect of magnetic buoyancy, et cetera, as was assumed by several astronomers, the inclb .ation of t11.e sun spot axis

y"nl be in the direction of the solar rotation in the lower latitude while this in

the opposite direction in the higher latitude. If we further imagine the

pos-sibility of twist of the horizontal magnetic lines of force by the effect of the shear, the twist wiU be most exaggerated near the region at about 40 degrees of latitude since the shear changes sign there and hence the rate of twist may become the largest. This situation may be interesting if we remember the

fact that these latitudes correspond roughly to the upper boundary of the Slm

spot belt. Similar curves are shcwl1 in Fig. 7 and 9 for the case with 5=-01>

and

R

t

=

0 18 7 S- .

The meridional current is shown in Fig. 10 which gives 'che

"- ,...

variation of U& versus latitude with

r,

as the parameter for the case with

S"" I.'S" and

Rl

=

01975'. FrOln these curves, we can conc1ude as follows:

The meridional current in the case wit11. radial viscosity larger tl~an the

(37)

the surface of the sun, falls along the equator and finally flows up to the higher latitude to complete a closed contour in the meridional plane. The trend of

other stream lines enclosed by the above peripheral stream line is similar .

The direction of the meridional current agrees with that expected from ob-servations. We do not find the type of meridional current distribution which

was inferred by Waldmeier from the sun spot motion (ReL 7) . . i~L·e. a dual meridional current distribution in which the low and the high latitude parts have opposite directions of circulation. Such a feature may be contained in the

higher order approximations of our theory.

Since the meridional current is proportional to

s-I ,

the casf'

with

s=

0#5'

may easily be obtained by only changing the nirection of the current.

Finally, since the rotational motion in the inner stable region

is a relic of the past related to the processes which oc.curred when the solar

system was first formed, our assumption of solid body rotation may be judged only af ter_ the establishment of the theory of the formation of the solar system. Sin1,~' ~lOwev~r, t~e process is exp~ct~d to ?e :rery violent, there may be room

for

'

&

hlghly diffuslve phenomenon wlthm a.lifetlme of the order of 10 7 to 108

yeaiis taking place in the early history of the sun. From such a point of view,

we ~lieve that our assumption does not contradict Cowling' s point of view.

~~. ;,.if :;:

(38)

1. 2. 3. 4. 5. 6. 7. 8. 9. Allen, C. W. Babcock, H. W. Biermann, L. Cowling, T. G. Kippenhahn, R.

"

Stromgren, B. Waldmeier, M. Ward, F. Wa siutyn ski, J. REFERENCES

Astro{:lhysical Quantities, 2nd Edition, University of London, The Athlone Press, 1963, p. p. 179

Astrophysical J., 133 (1961) p. p. 572

Zs. f. Ap., 28 (1951) p. p. 304

Solar Electrodynamics, in "The Sun" edited by Kuiper, G.P., The University of Chicago Press, (1953) p. p. 550-557 Astrophysical J., 137 (1963) p. p. 664 The

· .:sttn

.as.

a

~ar}n "The Sun" ibid p. p. 36-85

Ergebnisse und Problem der Sonnen forschung, 2nd Ed., (1955) p. p. 198 Astrophysical J., 141 (1965) p. p. 534

Studies in hydrodynamic s and structure

of stars and planets, Astrophysica Norvegica, Vol. 4, Oslo 1946 (to be obtained from Jacob Dybwad) p. p. 26-33

(39)

APPENDIX 1

The problem of solving equations like Eq. (118) and Eq. (119) satisfying the

boundary condition (53) can be generalized to the following problem: Solve

L(A

subject to the boundary condition --.;

t.f is analytic at

Equation (Al) is easily integrated to give:

=

d'r.

(m + 1)

(1S -

1) ;. 1

where the integration constant is determined to be zero by the boundary condition (A2). Finally, the integration of (A3) gives the following:

"V

U

=

_(m

_-+

_{1) [ (m ;. 1) ("{ - 1)} ₊_{1 }} + const. 27 (Al) (A2) (A3) (A4)

(40)

APPENDIX 2

~~II

=

B4 . F1 (1,

r; )

+ B5' F1 (2,

~

) - B6' F1 (3,

Ft ),

,...,

,...

U

eo , = Cl - C2 F2 (2,

r,

F2 (3,

r,

,..,

- C3

.

) + C4 . F2 (4,

r,

),

-

-~,

[

-

-P()II =

r,

C5'F1(1,r, ) + C6

.

F1 (2, (, ) '"

-

₎

}

+ C7

.

F1 (3,

r,

) + C8 F1 (4,

r,

-...,

-

,...

(j~,a. = U~I2.1 - ..1..

,

U~I/

-z..t

.L ""

l.{~IO

,

"

-'"

,...." Ucpl'2.l= Dl F1 (1,

r ,

) + D2 . F1 (2, "', "'- " -- D3 F1 (3,

r,

) - D4

.

F1 (4,

r,

) " -+ D5 F1 (5~

r,

),

" , U&02.1

=

D50 - D6 " - " -F2 (2,

r , ) -

D7 . F2 (3,

r,

)

- D8 ,.., " -F2 (4,

r, )

+ Dg . F2 (5,

n

-- D10 . F2 (6,

n ),

" - ,..,

-

-Wfj2.

=

'P02.1

-l~ol

_a

- .2:.

~eo 2'1 ~ _- _A!

[

--

. F1 (2,

r;

POf

I

=

r,

El F1 (1, t", ) + E2 ,....

-+ E3 . F1 (3,

r. )

+ E4 . F1 (4,

r,

)

+ E5 • F1 (5,

r; )

+ E6 . F1 (6,

r. )]

(A5) (A6) (A7) (A8) (Ag) (AlO) (All) (A12) (A13)

(41)

where

GM1 (A, B) = B

l

(B + 1) • A - B } ,

GM2 (A, B) = GM1 (A, B) . GM1 (A, B + 1),

F1 (A, B) F2 (A, B) Al =

Y

A2 = Al - 1, A3

=

~1 /.~ A4 =

Ja'

A2LAl, A5 =

J,.

A2/ (2.A1 - 1), A6 = A5' (1

+

A4), A7 =

j~.

A2, A8

=

A5

{J;.

A6 . A2 . (B12 - BU) -A9

=

A2· (1

+

A4) / 24, B11

=

1 GM1 (Al, l) , ;;..:: B12

=

1 , GM1 (Al, 2) B13

=

1 GM1 (Al, 3) B14

₌

1 GM1 (Al, 4) B15

=

1 -GM1 (Al, 5) B16

=

1 GM1 (Al, 6) 29

~

(1 +. A4)} (A14) (A15) (A16) (A17) (A18) (A19) (A20) (A21) (A22) (A23) (A24) (A25) (A26) (A27) (A28) (A29) (A30) (A31) (A32)

(42)

Cll

=

1 GM2 (Al, 1) , C12

=

1 GM2 (Al, 2) , C13

=

~.J,.

(A2.>2 . A6, D 11

=

J

_{7 ·} A2· A8, D12

=

J,.

J

₁2.· (A2)3. A6 , D13

=

D14 = D15

=

D11 . A2 J D12· A2,

J, .

(A2)2 _{A6/6 ,} B1

=

B7 . (B12 - B11) , B2

=

B8. B11 J B3 = B8· B12, B4

=

-[~.

A4 + B9 • (Bll - B12) } B5

=

B9· B 11 /2 B6

=

B9' B12/3 B7

=

Al A6 B8

=

A7 A6 B 9

=

J

₇ • A6 . A2 Cl

=

A8· Al· (B 12 - B 11) (A33) (A34) (A35) (A36) (A37) (A38) (A39) (A40) (A41) (A42) (A43) (A44) (A45) (A46) (A47) (A48) (A49)

+

C13· Al . {C11 (B13 - B11) - C12 . (B14 - B11)] (A50) C2

=

A8 C3

=

C13· C11 C4

=

C13· C12 C 5

=

A 7· B 11· [A8 + C 13· (C 11 - C 12) } (A51) (A52) (A53) (A54)

(43)

C 6

= -

A 7' A8 . B 12 • C7 = -A7· C13 • Cll • B13. , C8 = A7 . C13 . C12 . B14. Dl = Dll· (B12 - Bll)

+

D12 . [Cll' (B13 -B11) - C12 . (B14 - Bll)}

_•

D2 = [D11

+

D12 . (C11 - C12) } . Bll/2.0. D3 = Dll' B12/3. D4 = D12 . ·Cll· B13/4. D5 = D12· C12 . B14/5 •

d-2.

=

-Al [ A 9 . B11· (B12 - Bll) - D 13· [ B 12 [ B 13 . (B 14 - B 11) - B11 . (B12 - Bll

)J

- B11 .

l

B12. (B13 - B11) - B11. (B12 - B11

)}J

- D14. [C11.

l

B13· (B14. (B15 - B11) - B11· (B12 - B11» - B 11· (B 12. (B 13 - B 11 ) - B11· (B12 - Bll»)} - C12. { B14· (B15. (B16 - B11) - B11· (B12 - B11») - BIl· {B12. (B13 - B11) - B 11. (B 12 - B 11»)} ] 31 (A55) (A56) (A57) (A58) (A59) (A60) (A61) (A62)

(44)

" --D15 . [Bll . [B12. (B13 - B11) - Bll· (B13 - Bll)} - B12 . [ B13· (B14 - B11) - B 11· (B 12 - B 11) }

J ]

(A63) D50

=

-

J,

'te.

(A64) D6

=

~,.

B11 .[ Ag

+

D13· (B12 - B11)

+

D15· (B11 - B12)

+

D14 .

L

Cll· (B13 - B11) - C13· (B14 - Bll)E (A65) D7 =

J~.

B11· B12.[ D13 - D15

+

D14· (Cll - C12)

J

(A66) D8

=

J

_{ro •} B12· B13· (D15 - D13), (A67) Dg

=

J,.

Bl3· Bl4 • C11· D14, (A68) D10

=

J".

Bl4 . B15· C12 . Dl4 , (A69) El

= -

J,.

A7· B11 . [ - Ag . B11

+

D13· { B12· (B13 - B1l) - Bll . (B12 - B11)]

+

D14 .[C11 .[B13. (Bl4 - B11) - B11· (B12 - B11)] - C12 .[B14. (B15 - Bl1) - Bll· (B12 - Bll)}]

+

D15 . [B11. (B12 - B1l) - B12· (B13 - B1l)

1]

(A70) E2

=

-

J,.

A7 . B11 . B12 . [Ag

+

D13· (B12 - Bll)

+

D14 . [ Cll· (B13 - Bll) - C12 . (B14 - B11)} + D15· (Bll - Bl2)] E3

=

-J,.

A7 . Bll· B12 . Bl3

.

l

D13

+

D14 . (C11 - C12) - D15} (A7l) (A72)

(45)

E4 =

J,.

A7. B12 B13· B14. (D13 - D15) E5

=

JI .

A7. D14 Cll. B13 • B14. B15, EG

= -

J, .

A7. D14 • C12. B14. B15. B16 33 (A73) (A74) (A75)