[mesa-users] Rotation : relax_initial_omega_div_omega_crit + preserving omega in several inlists
Héctor MR
hector.mr at pitt.edu
Mon Jun 20 18:31:37 EDT 2016
Hello Frank,
Okay, it has been difficult to perform for a MESA inexpert like me but I
think have managed to do it! I would like to know if my explanation makes
sense. Sorry the email will be extremely long, but I have to include all
the code :) . It would be great if Bill could go over this; thanks in
advance to you guys!
First, I tried with extras_finish_step in run_star_extras.f but was
unable to write the parameters I wanted, as they do not belong to the
"star_info" structure. Then I realized I could not "change" MESA by adding
the write statements you suggested unless I executed
$MESA_DIR/star/text/mkx
and then $MESA_DIR/mk.
I did so in relax_omega_check_model under /star/private/relax.f90:
!
* kind_of_relax = 0 => target = new_omega ! kind_of_relax = 1 =>
target = new_omega_div_omega_crit ! kind_of_relax = 2 => target =
new_surface_rotation_v *
if (kind_of_relax == 0) then
new_omega = target_value
else if (kind_of_relax == 1) then
call set_surf_avg_rotation_info(s)
new_omega = target_value*s% omega_crit_avg_surf
else if (kind_of_relax == 2) then
new_omega = target_value*1d5/(s% photosphere_r*Rsun)
else
write(*,2) 'bad value for kind_of_relax', kind_of_relax
stop 'relax_omega_check_model'
end if
this_step_omega = frac*new_omega + (1 - frac)*starting_omega
if (s% model_number > starting_model_number + num_steps_to_use +
50 .or. &
abs(s% omega(1) - new_omega) < 1d-4*new_omega) then
write(*,2) 'final step: wanted-current, current, wanted', &
s% model_number, new_omega-s% omega(1), s% omega(1),
new_omega
relax_omega_check_model = terminate
s% termination_code = t_relax_finished_okay
else
write(*,2) 'relax to omega: wanted-current, current, wanted', &
s% model_number, new_omega-s% omega(1), s% omega(1),
new_omega
* write(*,'(A,3X,3(I3,1X))') 'ipar: starting_model_number,
num_steps_to_use, kind_of_relax', ipar
write(*,'(A,3X,2(ES16.6,2x))') 'rpar: target_value, starting_omega',
rpar write(*,'(A,1X,I2)') 'kind_of_relax',
kind_of_relax write(*,'(A,1X,ES16.6)') 'target_value',
target_value write(*,'(A,1X,ES16.6)') 'omega_crit_avg_surf', s%
omega_crit_avg_surf write(*,'(A,1X,ES16.6)') 'target_value *
omega_crit_avg_surf', target_value*s% omega_crit_avg_surf
write(*,'(A,1X,ES16.6)') 'new_omega', new_omega
write(*,'(A,1X,ES16.6)') 'photosphere_r', s% photosphere_r
write(*,'(A,1X,ES16.6)') 'target_value*1d5/(s% photosphere_r*Rsun)',
target_value*1d5/(s% photosphere_r*Rsun)*
end if
Then, I tried with both *set_initial_omega_crit*:
relax to omega: wanted-current, current, wanted 1
1.7908379232913815D-05 0.0000000000000000D+00 1.7908379232913815D-05
ipar: starting_model_number, num_steps_to_use, kind_of_relax
0 10 2
rpar: target_value, starting_omega 1.000000E-01
0.000000E+00
*kind_of_relax 2*
*target_value 1.000000E-01*
omega_crit_avg_surf 8.732251E-01
target_value * omega_crit_avg_surf 8.732251E-02
*new_omega 1.790838E-05*
photosphere_r 8.023188E-03
* target_value*1d5/(s% photosphere_r*Rsun) 1.790838E-05*
And *relax_surface_rotation_v*:
ipar: starting_model_number, num_steps_to_use, kind_of_relax
0 10 1
rpar: target_value, starting_omega 1.000000E-01
0.000000E+00
*kind_of_relax 1*
*target_value 1.000000E-01*
*omega_crit_avg_surf 8.732251E-01*
*target_value * omega_crit_avg_surf 8.732251E-02*
new_omega 8.732251E-02
photosphere_r 8.023188E-03
target_value*1d5/(s% photosphere_r*Rsun) 1.790838E-05
So that both instructions seem to be interchanged due to kind_of_relax. I
checked /star/run_star_support.f90 and:
if (*s% rotation_flag .and. s% job% relax_omega*) then
write(*,1) 'new_omega', s% job% new_omega
*call star_relax_uniform_omega( &*
*id, 0*, s% job% new_omega, s% job%
num_steps_to_relax_rotation,&
s% job% relax_omega_max_yrs_dt, ierr)
if (failed('star_relax_uniform_omega',ierr)) return
end if
if (*s% rotation_flag .and. s% job% relax_initial_omega .and.
.not. restart*) then
*call star_relax_uniform_omega( &*
*id, 0*, s% job% new_omega, s% job%
num_steps_to_relax_rotation,&
s% job% relax_omega_max_yrs_dt, ierr)
if (failed('star_relax_uniform_omega',ierr)) return
write(*,1) 'new_omega', s% job% new_omega
end if
if (*s% rotation_flag .and. s% job% relax_surface_rotation_v*) then
*call star_relax_uniform_omega( &*
*id, 1*, s% job% new_surface_rotation_v, s% job%
num_steps_to_relax_rotation,&
s% job% relax_omega_max_yrs_dt, ierr)
if (failed('star_relax_uniform_omega',ierr)) return
s% job% new_omega = s% job% new_surface_rotation_v*1d5/s% r(1)
write(*,1) 'new_surface_rotation_v', &
s% job% new_surface_rotation_v, s% job% new_omega
end if
if (
*s% rotation_flag .and. & s% job%
relax_initial_surface_rotation_v .and. .not. restart*) then
write(*,1) 'new_omega', s% job% new_omega
write(*,*) 'call star_relax_uniform_omega'
*call star_relax_uniform_omega( &*
*id, 1*, s% job% new_surface_rotation_v, s% job%
num_steps_to_relax_rotation,&
s% job% relax_omega_max_yrs_dt, ierr)
if (failed('star_relax_uniform_omega',ierr)) return
write(*,2) 'new_surface_rotation_v', &
s% model_number, s% job% new_surface_rotation_v
end if
if (*s% rotation_flag .and. s% job% relax_omega_div_omega_crit*)
then
if (failed('star_surface_omega_crit',ierr)) return
*call star_relax_uniform_omega( &*
*id, 2*, s% job% new_omega_div_omega_crit, &
s% job% num_steps_to_relax_rotation,&
s% job% relax_omega_max_yrs_dt, ierr)
if (failed('star_relax_uniform_omega',ierr)) return
write(*,2) 'new_omega_div_omega_crit', &
s% model_number, s% job% new_omega_div_omega_crit
end if
if (
*s% rotation_flag .and. & s% job%
relax_initial_omega_div_omega_crit .and. .not. restart*) then
if (failed('star_surface_omega_crit',ierr)) return
*call star_relax_uniform_omega( &*
*id, 2*, s% job% new_omega_div_omega_crit, &
s% job% num_steps_to_relax_rotation,&
s% job% relax_omega_max_yrs_dt, ierr)
if (failed('star_relax_uniform_omega',ierr)) return
write(*,2) 'new_omega_div_omega_crit', &
s% model_number, s% job% new_omega_div_omega_crit
end if
With the corresponding subroutine under /star/public/star_lib.f90:
*subroutine star_relax_uniform_omega(id, &*
*kind_of_relax*, target_value, num_steps_to_relax_rotation, &
relax_omega_max_yrs_dt, ierr)
use relax, only: do_relax_uniform_omega
integer, intent(in) :: id, kind_of_relax,
num_steps_to_relax_rotation
real(dp), intent(in) :: target_value,relax_omega_max_yrs_dt
integer, intent(out) :: ierr
call do_relax_uniform_omega(id, &
kind_of_relax, target_value, num_steps_to_relax_rotation, &
relax_omega_max_yrs_dt, ierr)
end subroutine star_relax_uniform_omega
What is your opinion on this matter? Am I raving or is it correct?
Thanks!
--
Héctor
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.mesastar.org/pipermail/mesa-users/attachments/20160620/88c54b23/attachment.html>
More information about the Mesa-users
mailing list