2023.08.31

第19回 HPC-Phys 勉強会

ラグランジュ的定式化を⽤いた
2次元の星の平衡形状の新たな計算法
早稲⽥⼤学 ⼩形 美沙

Contents
p Introduction
p Numerical Methods
ü Conﬁguration (Lagrangian formulation)
ü Self-gravity (Spectral method)p Results
p SummaryIntroduction
n 1D stellar evolution calculation
Basic calculation
to understand stellar evolution15)

105log(L/(Henyey et al. 1959)203

ü Angle dependence averaged11Ω/Ωcrit=0.0
0.1
0.3
0.5
0.7
0.8
0.9
0.99
1.0

0.9
0.8
0.7

z0.6log(T / )0.5
0.4
0.30.2

Roche model

0.1
00

0.2

0.4

0.6

0.8

1

1.2

1.4x

2023.08.31

第19回 HPC-Phys 勉強会

Introduction
Calculations correctly evaluate the eﬀect of rotation• All previous stellar evolution calculations are 1D.
• No 2D structure calculations using Lagrangian
coordinate.

2023.08.31

第19回 HPC-Phys 勉強会

Introductionn Eulerian/Lagrangian coordinates• Eulerian︓Coordinates fixed in space

Solve for the value of a physical quantity on coordinate• Lagrangian︓Coordinates fixed to a fluid element
Solve for where the physical quantity moves to2D Eulerian coordinate

2D Lagrangian coordinate

Solve for the value of a physical quantity

Solve for the position
2023.08.31

第19回 HPC-Phys 勉強会

Introduction
ü Suitable coordinates for the expansion/contraction of stars
ü Appropriate introduction of rotation eﬀects➡ 2D Lagrangian coordinate

1D Lagrangian
coordinate

2D Euleriancoordinate

2D Lagrangian
coordinate
2023.08.31

第19回 HPC-Phys 勉強会

Introductionn Flow of evolution calculation
Dynamicaltimescale≪

Evolution
timescale

Ⅰ. Hydrostatic equilibrium calculation
Ⅱ. Time evolution calculation
Nuclear network
Energy transportation
Repeat ➡ Evolution Calculation
2023.08.31

第19回 HPC-Phys 勉強会

Contents
p Introduction
p Numerical Methods
ü Configuration (Lagrangian formulation)
ü Self-gravity (Spectral method)p Results
p SummaryFlow of calculation
ConfigurationüLagrangian formulation
üW4 method
Self gravity
ü Spectral method
Remeshingü Interpola7on
2023.08.31

第19回 HPC-Phys 勉強会

Method: ConfigurationnLagrangian formulation
Force-balance eq.︓
1
1 !
𝛁P = −𝛁𝜙 + Ω 𝛁 𝑅 sin Θ !
𝜌2Solve for the position(Not for values of physical quantities)↓If not set up appropriately,
no solution exists.

Assumption: axisymmetry, equatorial symmetry
2023.08.31

第19回 HPC-Phys 勉強会

Method: Configurationn Flow of conﬁguration calculation
1. Setting distributions of mass, entropy and angularmomentum
2. Finding the hydrostatic equilibrium conﬁguration withoutchanging the given physical quantities
Balance eq.︓

1
1 !
𝛁P = −𝛁𝜙 + Ω 𝛁 𝑅 sin Θ !
𝜌2mass, entropy (pressure)angular momentum

2023.08.31

第19回 HPC-Phys 勉強会

Method: Configurationn Two types of coordinates
ü Lagrangian coordinate
Spherical reference configuration＝Give the distribution of physical quantities

ü Eulerian coordinate
Non-spherical configuration＝Spatial position

2023.08.31

第19回 HPC-Phys 勉強会

Method: Configurationn Flow of conﬁguration calculation
1. Setting distributions of mass, entropy and angularmomentum
2. Finding the hydrostatic equilibrium configuration withoutchanging the given physical quantities
Balance eq.︓
mass, entropy (pressure)1
1 !
𝛁P = −𝛁𝜙 + Ω 𝛁 𝑅 sin Θ !
𝜌2angular momentum

EOS: 𝑃 = 𝑃 𝜌, 𝑠Rotation law: Ω = Ω 𝑅, Θ*+*,Δ𝑉 = det 𝐽 Δ𝑣Density: 𝜌 =𝐽=

𝜕𝑅 !/3#𝜕𝑟 !/3

𝜕 cos Θ)𝜕𝑟 !/3 "#$ %𝜕𝑅 !/3#
𝜕 cos 𝜃 !

𝜕 cos Θ)𝜕 cos 𝜃 &!/!

2023.08.31

第19回 HPC-Phys 勉強会"#$ %& /!

Method: Configuration
Lagrangian coordinate(Spherical reference conﬁgura5on)
Eulerian coordinate

11 !
𝛁P + 𝛁𝜙 − Ω 𝛁 𝑅 sin Θ ! = 0 ︓Balance eq.
𝜌2
𝜕𝑃
𝜕𝜙+𝜌
− 𝜌 𝑅 sin Θ Ω" sin Θ = 0
𝜕𝑅𝜕𝑅𝜕𝑃
𝜕𝜙𝐹# =+𝜌− 𝜌 𝑅 sin Θ Ω" cos Θ = 0𝑅𝜕Θ𝑅𝜕Θ𝐹! =𝐹&=𝐹%𝜕𝑅𝜕𝑟
𝜕𝑅𝑟𝜕𝜃

R𝜕Θ𝜕𝑟
𝑅𝜕Θ𝑟𝜕𝜃

𝐹(𝐹)︓coordinate 𝑅, Θ → 𝑟, 𝜃𝑑𝑃
𝜕𝑅𝑅𝜕Θ𝜌𝑗 "𝜕𝑅𝑅𝜕Θ!
#𝐹$ =+𝜌
𝐹%&'( +𝐹%&'( +sinΘ
+cos Θ = 0𝑑𝑟𝜕𝑟
𝜕𝑟
𝑅 sin Θ ) 𝜕𝑟
𝜕𝑟𝑑𝑃
𝜕𝑅𝑅𝜕Θ
𝜌𝑗 "
𝜕𝑅𝑅𝜕Θ
!
#𝐹* =+𝜌
𝐹%&'( +𝐹%&'( +sinΘ
+cos Θ = 0𝑟𝑑𝜃𝑟𝜕𝜃
𝑟𝜕𝜃
𝑅 sin Θ ) 𝑟𝜕𝜃
𝑟𝜕𝜃
2023.08.31

第19回 HPC-Phys 勉強会

Method: ConfigurationDiﬀerentiation︓
𝐹" =𝑃#$ + 𝑃## − 𝑃$$ − 𝑃$#+𝜌𝑟% − 𝑟$𝑅%# − 𝑅$#𝑅## Θ%# − Θ$#*𝐹&'() +𝑟% − 𝑟$
𝑟% − 𝑟$+𝐹, =𝑃$# + 𝑃## − 𝑃$$ − 𝑃#$+𝜌
𝑟# 𝜃% − 𝜃$𝜌𝑗 #𝑅## sin Θ## %𝑅#% − 𝑅#$𝑟# 𝜃% − 𝜃$+

𝑅#% − 𝑅#$𝑅## Θ%# − Θ$#sin Θ## +𝑟% − 𝑟$
𝑟% − 𝑟$*𝐹&'()+

𝜌𝑗 #𝑅## sin Θ## %+𝐹&'()𝑅## Θ#% − Θ#$𝑟# 𝜃% − 𝜃$

𝑅#% − 𝑅#$𝑟# 𝜃% − 𝜃$cos Θ## = 0+𝐹&'()sin Θ## +𝑅## Θ#% − Θ#$𝑟# 𝜃% − 𝜃$𝜌=cos Θ## = 0

𝜌** + 𝜌*+ + 𝜌+* + 𝜌++
4fixed2023.08.31

第19回 HPC-Phys 勉強会

Method: Configurationn W4 method

Okawa et al. 2023

Root-ﬁnding scheme for nonlinear simultaneous equations
using iterative methods• Local convergence similar to the Newton method• Excellent global convergence→ Useful for stiﬀ equations which cannot be solvedby the Newton method𝑑" 𝒙𝑑𝒙+𝑀+ 𝑀" 𝑭 𝒙 = 0+𝑑𝜏 "𝑑𝜏⟹4𝑑𝒙= 𝑋𝒑𝑑𝜏𝑑𝒑= −2𝒑 − 𝑌𝑭𝑑𝜏𝒙"#$ = 𝒙" + Δ𝜏𝑋𝒑"𝒑"#$ = 1 − 2Δ𝜏 𝒑" − Δ𝜏𝑌𝑭 𝒙"
2023.08.31

第19回 HPC-Phys 勉強会

Flow of calculation
ConfigurationüLagrangian formulation
üW4 method
Self gravity
ü Spectral method
Remeshingü Interpola7on
2023.08.31

第19回 HPC-Phys 勉強会

Method: Self-gravity
l Solve Poisson's equation by spectral method with Legendre
polynomial 𝜇 = cos ΘUse the Eulerian coordinateü Laplacian takes the simplest form
ü Application of the spectral method becomes easiest∇" 𝜙 𝑅, 𝜇 = 4𝜋𝐺𝜌 𝑅, 𝜇 →Eulerian coordinate

𝜕
𝜕𝜙𝜕
𝜕𝜙𝑅"
+1 − 𝜇"= 4𝜋𝐺𝜌𝑅 "𝜕𝑅𝜕𝑅𝜕𝜇𝜕𝜇,!"# /!"#

𝜙 𝑅, 𝜇 = R R 𝒂𝒍𝒎 𝑃T, 𝑅 𝑃/ 𝜇,-. /-.
,!"# /!"#4𝜋𝐺𝜌 𝑅, 𝜇 𝑅 " = R R 𝒃𝒍𝒎 𝑃T, 𝑅 𝑃/ 𝜇,-. /-.

𝑃/ : Legendre polynomial
𝑃T, : shifted Legendre polynomial𝑅∗2023.08.31

第19回 HPC-Phys 勉強会

Method: Self-gravity
l Solve Poisson's equation by spectral method with Legendrepolynomial∇" 𝜙 𝑅, 𝜇 = 4𝜋𝐺𝜌 𝑅, 𝜇 →

𝜕𝜕𝜙𝜕
𝜕𝜙𝑅"
+1 − 𝜇"= 4𝜋𝐺𝜌𝑅 "𝜕𝑅𝜕𝑅𝜕𝜇𝜕𝜇Legendre polynomial𝑑𝑑
1 − 𝑥"
𝑓 𝑥𝑑𝑥𝑑𝑥+𝜆 𝜆+1 𝑓 𝑥 =0

Shifted Legendre polynomial

𝑥+12!∗𝑅∗ZZY 𝑃, 𝑅 𝑃/ 𝑅 𝑑𝑅 =𝛿,/2𝑙+1.
𝑑2 𝑑Z
𝑃 𝑅 =𝑃 𝑥𝑑𝑅 ,Δ𝑅 𝑑𝑥 ,−1,1 → 0, 𝑅∗

𝑅 = 𝑅∗

2023.08.31

第19回 HPC-Phys 勉強会

Method: Self-gravity
𝜕
𝜕𝜙𝜕
𝜕𝜙
""
𝑅+1−𝜇= 4𝜋𝐺𝜌 𝑅, 𝜇 𝑅 "𝜕𝑅𝜕𝑅𝜕𝜇𝜕𝜇,!"# /!"#

𝜙 𝑅, 𝜇 = R R 𝑎,/ 𝑃T, 𝑅 𝑃/ 𝜇,-. /-.
,!"# /!"#

4𝜋𝐺𝜌 𝑅, 𝜇 𝑅 " = R R 𝑏,/ 𝑃T, 𝑅 𝑃/ 𝜇,-. /-.

RR𝑑𝑑
𝑅"𝑃T, 𝑅𝑑𝑅𝑑𝑅𝑎,/ 𝑃/ 𝜇 − 𝑚 𝑚 + 1 𝑎,/ 𝑃T, 𝑅 𝑃/ 𝜇𝑑𝑑"R𝑅𝑃T, 𝑅𝑑𝑅𝑑𝑅𝑎,/ − 𝑚 𝑚 + 1 𝑎,/ 𝑃T, 𝑅= R 𝑏,/ 𝑃T, 𝑅 𝑃/ 𝜇= R 𝑏,/ 𝑃T, 𝑅

2023.08.31

第19回 HPC-Phys 勉強会

Method: Self-gravity
n Diﬀerentiation of Legendre polynomial
𝒅2 𝑑T
𝑷𝒍 𝑹 =𝑃, 𝑥 = R 𝐶,,% 𝑃Z ,% 𝑅𝒅𝑹
Δ𝑅 𝑑𝑥𝑑𝑃 𝑥 = 2𝑙 − 1 𝑃,3+ 𝑥 + 2𝑙 − 5 𝑃,3) 𝑥 + ⋯𝑑𝑥 ,0102𝐶%%, =Δ𝑅030105030⋱0105…
…
…

⋮02𝑙 − 50
2𝑙 − 102023.08.31

第19回 HPC-Phys 勉強会

Method: Self-gravity
n Multiplying Legendre polynomial by 𝑅 F(Arfken 1985, p. 700)g
𝑅 " 𝑃T, 𝑅 = R 𝐷,4 𝑃
4 𝑅+++2𝑛 + 1 Δ𝑅 )𝑎"
𝐷,4 =Y 𝑃, 𝑥 𝑃4 𝑥 𝑑𝑥 + 2𝑎 Y 𝑥𝑃, 𝑥 𝑃4 𝑥 𝑑𝑥 + Y 𝑥 " 𝑃, 𝑥 𝑃4 𝑥 𝑑𝑥
Δ𝑅
22𝑛 + 1 3+3+3+2𝑛 + 1 Δ𝑅 )𝑎"+"=𝛿,4 + 2𝑎𝐷,4 + 𝐷,4Δ𝑅
22𝑛 + 1*

𝐷,4+ :"
𝐷,4 :4 𝑥𝑃. 𝑥 𝑃/ 𝑥 𝑑𝑥 =-**4 𝑥 +𝑃. 𝑥 𝑃/ 𝑥 𝑑𝑥 =-*𝑎=

𝑅56 + 𝑅789
Δ𝑅2 𝑙+1for 𝑛 = 𝑙 + 12𝑙 + 1 2𝑙 + 32𝑙
for 𝑛 = 𝑙 − 12𝑙 − 1 (2𝑙 + 1)

2 𝑙+1 𝑙+22𝑙 + 1 2𝑙 + 3 2𝑙 + 52 2𝑙+ + 2𝑙 − 12𝑙 − 1 2𝑙 + 1 2𝑙 + 32𝑙 𝑙 − 12𝑙 − 3 2𝑙 − 1 2𝑙 + 1for 𝑛 = 𝑙 + 2for 𝑛 = 𝑙
for 𝑛 = 𝑙 − 22023.08.31

第19回 HPC-Phys 勉強会

Method: Self-gravity
R𝜕𝜕
𝑅"𝑃T, 𝑅𝜕𝑅
𝜕𝑅𝑎,/ − 𝑚 𝑚 + 1 𝑎,/ 𝑃T, 𝑅= R 𝑏,/ 𝑃T, 𝑅

R R R R 𝐶,,%%% 𝐷,%%%,%% 𝐶,%%,% − 𝑚 𝑚 + 1 𝛿,,% 𝑎,%/ 𝑃T, 𝑅 = R 𝑏,/ 𝑃T, 𝑅= 𝑴j𝒍𝒍2
R R 𝑀,,: % 𝑎,%/% 𝑃T, 𝑅 = R 𝑏,/% 𝑃T, 𝑅 → R 𝑀,,: % 𝑎,%/% = 𝑏,/%𝑎,%/% = 𝑀,,3+% 𝑏,/%𝑴j𝒍𝒍2 =0𝑚& -line

02023.08.31

第19回 HPC-Phys 勉強会

Method: Self-gravity
n Boundary conditions
ü Inner boundary
𝜕𝜙/%k=0𝜕𝑅 !-.

𝜙/% = R 𝑎,/% 𝑃Z, 𝑅,⟹ R 𝑎,/% × −1 , 𝑙 𝑙 + 1 = 0,

ü Outer boundary𝜕𝜙/%𝑚: + 1k=−𝜙/%𝜕𝑅 !-!
𝑅∗∗

𝑴j𝒍𝒍2 =00

⟹ R 𝑎,/% × 𝑚: + 1 + 𝑙 − 1 𝑙 = 0,𝑚& -line:
Outer boundary𝑚'() -line:
Inner boundary
2023.08.31

第19回 HPC-Phys 勉強会

Method: Self-gravity
n Coeﬃcient of density.-./ 2-./4𝜋𝐺𝜌 𝑅, 𝜇 𝑅 + = F F 𝒃𝒍𝒎𝑃H. 𝑅 𝑃2 𝜇.01 201

ü Gaussian quadrature method
2𝑙 + 1 2𝑚 + 1 !∗ +𝑏,/ =Y Y 4𝜋𝐺𝜌 𝑅, 𝜇 𝑅 " 𝑃Z, 𝑅 𝑃/ 𝜇 𝑑𝑅𝑑𝜇𝑅∗2.
3+4&𝑃/ 𝑦;' 𝑤;& 𝑤;';& -+ ;'-+𝜌"+𝜌++

4'2𝑙 + 1 2𝑚 + 1=
R R 4𝜋𝐺𝜌 𝑥;& , 𝑦;' 𝑥;"& 𝑃T, 𝑅 𝑥;&𝑅∗2𝜌""𝜌+"))s< 𝛼 𝑀s= 𝛽 𝜌 𝑟<= , 𝜃<=𝜌: 𝑟′, 𝜃′ = R R 𝑀<-+ =-+

I* 𝛼 = − 5 1 − 𝛼
𝑀+L
𝑀+ 𝛼 = 1 + 𝛼 1 − 𝛼
I! 𝛼 = 5 1 + 𝛼
𝑀−1 ≤ 𝛼, 𝛽 ≤ 1+2023.08.31

第19回 HPC-Phys 勉強会

Method: Self-gravity
n Expansion coeﬃcients of gravitational force
𝑎%, *, = 𝑀%%+$, 𝑏%*,!
𝐹%&'(𝑅, 𝜇𝜕𝜙== R R R 𝐶,,% 𝑎,%/ 𝑃Z, 𝑅 𝑃/ 𝜇𝜕𝑅,!"# 3+ /!"#= R,-.

#𝐹%&'( 𝑅, 𝜇 =𝑹

R 𝑮𝒍𝒎 𝑃T, 𝑅 𝑃/ 𝜇/-.

𝜕𝜙sin Θ 𝜕𝜙== R R R 𝑎,/% 𝐶/%/ 𝑃Z, 𝑅 𝑃/ 𝜇𝑅𝜕Θ𝑅 𝜕 cos Θ,!"# /!"# 3+= R,-.

𝚯R 𝑮𝒍𝒎 𝑃T, 𝑅 𝑃/ 𝜇/-.

2023.08.31

第19回 HPC-Phys 勉強会

Flow of calculation
ConﬁgurationüLagrangian formula7onüW4 method
Self gravity
ü Spectral method
Remeshingü Interpola7on
2023.08.31

第19回 HPC-Phys 勉強会

Method: Remeshingn Remeshing

Dealing with problems due to mesh deformation(due to numerical errors associated with finite differences)
1. Create a new smooth grid
2. Redistribute conserved quantities on the new grid10.80.6z

ü Inaccuracy
ü Non-convergence

Zigzag

0.4

0.2

00

0.2

0.4

0.6

0.8

1

1.2

1.4

x

2023.08.31

第19回 HPC-Phys 勉強会

Contents
p Introduction
p Numerical Methods
ü Configuration (Lagrangian formulation)
ü Self-gravity (Spectral method)p Results
p SummaryModell EOS: Barotropeü Isobaric and isopycnic surfaces coincide$$#𝑃 = 𝐾, 𝜌𝑁 = 1.5l Rotation law
For 𝑟 sin 𝜃 ≫ 1, 𝑗-constant law (Eriguchi & Mueller 1985)𝑟, 𝜃𝑅, Θ𝑗. 𝑟 sin 𝜃 "𝑗 𝑟, 𝜃 =1 + 𝑟 sin 𝜃 "2023.08.31

第19回 HPC-Phys 勉強会

Results (Barotrope)n Models with diﬀerent rotation rates
Intermediate𝑅.
= 0.743
𝑅/0

1.6
Shrinkρ^

z1𝑊 ≡ ∫ 𝜌𝜙 𝑑𝑉
2

𝑅.
= 0.380
𝑅/0𝑇= 0.168
𝑊

𝑇= 5.46×10+!
𝑊1
0.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0

Gravitational energy :Rapid

Rotational energy :𝑇≡1∫ 𝜌Ω+ 𝑅 sin Θ + 𝑑𝑉
2

spherical
intermediate rotation
rapid rotation

1.4
1.2
1

Elongate

0.8
00.6

0.5

0.4 spherical

1

1.5

2

2.5

x2023.08.31
rapid rotation

第19回 HPC-Phys 勉強会

Results (Barotrope)n Angular momentum distribution
ü Angular momentum is determined solely by distancefrom the axis of rotation
𝑁 = 1.5
𝑁1 , 𝑁2 = 16, 16

Intermediate rotation model

1

0.6

0.9
0.5

0.8
0.7

0.4
0.3

j0.5

^

z

0.6

0.4
0.2

0.3
0.2

0.1

0.1
0000.20.4

0.6

0.8x11.2

1.42023.08.31

第19回 HPC-Phys 勉強会

Results (Barotrope)
Gravitational energy :n Convergence of accuracy
Our finite-difference scheme is of secondorder accuracy.
𝟑𝑼 + 𝑾 + 𝟐𝑻
𝑽𝒄 =𝑼 + 𝑾 + 𝑻-1

101𝑊 ≡ ∫ 𝜌𝜙 𝑑𝑉
2Rotational energy :𝑇≡1∫ 𝜌Ω+ 𝑅 sin Θ + 𝑑𝑉
2

Integrated pressure :𝑈 ≡ ∫ 𝑃 𝑑𝑉

𝑉3 ∝ 𝑁1+!

-2

Vc

10

10-3

N=1.0 (spherical)( rotation)
N=1.5 (spherical)( rotation)
N=2.5 (spherical)( rotation)

10-4

10-5

1020
Nr (Nt)30

402023.08.31

第19回 HPC-Phys 勉強会

Summary
A new method for calculating the hydrostatic equilibrium structureof stars in 2D using a Lagrangian formulation
Conﬁguration
ü Lagrangian formulation
ü W4 method
Self-gravity
ü Spectral method
Remeshingü Interpolation

Future Workü Multilayered stars
ü Introduction to various physics
ü Application to general relativity(Okawa et al. 2022)
2023.08.31

第19回 HPC-Phys 勉強会