# Pulsating stars: theory

We will review the classical models of star evolution and show how they relate to observations of variable stars.


## Historical developments (1900-1960)

Schwarzchild, Shapley, Eddington, Rosseland, Ledoux and Erika Vitense-Bohm are a few among the key scientists that helped shaping the theory of self-maintained oscillations in stars.

<a href="https://en.wikipedia.org/wiki/Paul_Ledoux"><img src="https://upload.wikimedia.org/wikipedia/commons/b/be/Paul_Ledoux.jpg"></a>

The paper published in the Zeitschrift fur Astrophysik by N. Baker and Kippenhahn in August 1961 was the crowning achievement of 60 years of intense research. These authors were the first to introduce the term of 'kappa mechanism' and to show the consistency of a complete physical model of stellar oscillations by numerical simulations. 

Rudolf Kippenhahn:

<a href="https://d.wattpad.com/story_parts/703530646/images/158c9a79dd2ab92b135113103519.jpg"><img src="https://d.wattpad.com/story_parts/703530646/images/158c9a79dd2ab92b135113103519.jpg"></a>


## Simplified model

We now will review an elementary model of pulsating star that produces a prediction of the fondamental frequency (or period) that we should expect.
We closely follow the lecture of Rolf Kudritzki (SS 2015).

Suppose that the density in the star is given by a function $\rho_0(r)$. Call $M_0(r)$ the mass of the star contained in the sphere
of radius $r$.

Assume that the star is in hydrostatic equilibrium. Then


\begin{equation} 
\begin{split}
\frac{dP_0}{d r_0} & = G M(r_0) \frac{\rho_0(r_0)}{r_0^2} \\
dM(r_0) & = 4  \pi  \rho_0(r_0) r_0^2  d_r0 \\
M(r_0) & = \int_0^{r_0} dM(r_0)
\end{split}
(\#eq:sim1)
\end{equation} 

The symbol $r_0$ refers to the position (radius) of material particules at equilibrium.

Let us now perturb this equilibrium and compress the star.

We assume that $r_0$ moves to  $r(t)$. We adopt a linear strain model in which strain is uniform in the star (hence strain does not depend on the radius).
The compression/deformation is homogeneous in the model:

\begin{equation} 
\begin{split}
r(t) &= r_0 (1+ x(t)) 
\end{split}
(\#eq:sim2)
\end{equation}

When we follow any given set of particles, the mass of this lagrangian set is conserved during deformation, therefore 

\begin{equation} 
\begin{split}
\rho r^2 dr & = \rho_0 r_0^2 dr_0 \\
r(t) &= r_0 (1+ x(t)) \\
dr &= (1+x)dr_0 \\
\rho = \rho_0(1+x)^{-3}
\end{split}
(\#eq:sim3)
\end{equation}


# Experimental confirmation


```{r}
au <- read.csv("https://mathstat.dal.ca/~fullsack/DATA/Auto.csv",
   header=T,quote="",sep=',')
au
```


```{r}
pl <- read.csv("https://mathstat.dal.ca/~fullsack/DATA/cepheidPL.csv",
   header=T,quote="",sep=',')
pl


x=pl$logPdays
y=pl$MinusMagV
plot(x,y,xlab="logPdays",ylab="-MagV",main="Relation Period-magnitude in Cepheids")
```



We would like to fit a straight line to these points.  We assume
a simple linear regression model
$$y_i = \beta_0 + \beta_1 x_i + \epsilon_i$$

and estimate the intercept $\beta_0$ and slope $\beta_1$, using
the method of least squares.  The least squares estimates $\hat \beta_0$ and
$\hat \beta_1$ are used to get predicted values for the outcome variable $y$,
as

$$\hat y_i = \hat \beta_0 + \hat \beta_1 x_i$$

The difference between the observed and predicted $y$'s are called residuals.
The $i$'th residual is

$$e_i = y_i - \hat y_i$$

and the sum of squared residuals is called the error sum of squares.

$$SSE = \sum_{i=1}^n (y_i - \hat y_i)^2$$

In R the calculations are as follows:

```{r}
lm.out=lm(y~x)
coef(lm.out) #least squares estimates
ypred=predict(lm.out,x=x) #predicted values
resids=y-ypred  #residuals
SSE=sum(resids^2) #error some of squares
print(paste("Error sum of squares is ",SSE))
```


```{r}
plot(x,y, main="Data plot with least squares line added")
abline(lm.out)
```


```{r}
plot(x,y, xlim=c(0,2),ylim=c(0,6),main="Data plot with least squares line added")
abline(lm.out)
```


Is the intercept null?


```{r}
summary(lm.out)

```

This suggests that the slope is of order 3.11


The simple model above suggests that the theoretical slope should be 3.3. Not bad!




