c  Mahc_Init_Data: Various Initial Data Routines for Mahc_Evolve
c  Copyright (C) 2001 Mark Miller

/*@@
  @file      tov_series.F
  @date      August 2000
  @author    Mark Miller
  @desc
     Calculates a series of TOV solutions.
  @enddesc
@@*/

#include "cctk.h"
#include "cctk_Parameters.h"
#include "cctk_Arguments.h"
#include "cctk_DefineThorn.h"

/*@@
  @routine    tov_series
  @date       August 2000
  @author     Mark Miller
  @desc
    Calculates a series of TOV solutions.
  @enddesc
  @calls
  @calledby
  @history
  @endhistory
@@*/


      subroutine tov_series(CCTK_FARGUMENTS)
      implicit none
      DECLARE_CCTK_FARGUMENTS
      DECLARE_CCTK_PARAMETERS

      CCTK_REAL ln_start, ln_end,tov_centden,rmaxtov, tov_eos_gamma
      CCTK_REAL tov_eos_k,star_mass,star_rad
      integer i,num,ntov,imaxtov
      write(*,*) '*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/'
      write(*,*) '*/*/   TOV_SERIES:'
      write(*,*) '*/*/      eos_k: ',eos_k
      write(*,*) '*/*/      eos_gamma: ',eos_gamma
      write(*,*) '*/*/      mass(m_dot) vs. centden dumped to fort.72 '
      write(*,*) '*/*/      radius vs. centden dumped to fort.73 '
      write(*,*) '*/*/      (radius)/(schwarz_radius) vs. centden'
      write(*,*) '*/*/                                to fort.74 '
      write(*,*) '*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/'
      ln_start = log(tov_series_start_centden)
      ln_end   = log(tov_series_end_centden)
      num = tov_series_num
      do i=0,num
c        constant log-centden spacing
         tov_centden = ln_start + (ln_end-ln_start)*dble(i)/dble(num)
         tov_centden = exp(tov_centden)
         rmaxtov = mahc_init_rmaxtov
         ntov = mahc_init_ntov
         tov_eos_gamma = eos_gamma
         tov_eos_k = eos_k
         call solvetov_iso(presstov,rhotov,epstov,mtov,phitov,r_schwartov,
     .      r_isotov,ntov,imaxtov,tov_centden,tov_eos_gamma,tov_eos_k,rmaxtov,
     .      mahc_init_eos_handle)
         star_mass = mtov(imaxtov)
         star_rad = r_isotov(imaxtov)
         write(72,*) tov_centden,'   ',star_mass
         write(73,*) tov_centden,'   ',star_rad
         write(74,*) tov_centden,'   ',star_rad/(2.0d0*star_mass)
      enddo

      stop

      return
      end