[mesa-users] solar model calibration

WU T. ggwutao at 126.com
Mon Aug 22 07:38:32 EDT 2016


Hi, everyone.
     I'm making a standard solar model via the 'pulse' package of v6208. In the process of generating it, we use the method of iteration to find the best fitting input parameters, mixing length parameter alpha_MLT, initial_Y, and initial_Z, when the star arriving the age of current solar age (4.57Gyr), the solar radius, luminosity, and the surface abundance ratio (Z/X)s match with observations 1.0+/-1e-5, 1.0+/-1e-5, and 0.0229+/-1e-5, respectively. Finally, I obtained the input parameters and the corresponding output parameters as follows:
The input parameters:
    alpha_MLT                                  initial_Y                             initial_Z
   1.827359206341649E+00    2.740249714286597E-01    1.894746401823611E-02
The final output parameters:
 L                                                  R                                       (Z/X)s                                         Rcz
1.000007508735689E+00    1.000004947388044E+00    2.289945347234215E-02    7.150509481272301E-01   
age                                              H1                                          H2                                           He3                               He4
4.570000000000000E+09  7.392063449196695E-01    2.432632994723156E-17    2.688521006888449E-05    2.438393485683135E-01

Based on the above initial inputs, we obtained a larger Rcz=0.71505 (the position of convective envelope), i.e. the convective envelope is shallower than observation. It is much larger than 0.713 (Basu & Antia 1997).

So, who can give some advices for this problem?
Thanks in advance!
cheers!
Tao Wu

inlist file
&star_job  ! HD49385  
      show_log_description_at_start = .false.     
      astero_just_call_my_extras_check_model = .true.
      create_pre_main_sequence_model = .true.

      pre_ms_relax_num_steps = 50     
      new_surface_rotation_v =0 ! solar (km sec^1)
      set_near_zams_surface_rotation_v_steps = -1 ! to turn on rotation when near zams
         ! if rotation_flag is false and L_nuc_burn_total >= L_phot and this control > 0
         ! then turn on rotation and set set_surf_rotation_v_step_limit to
         ! the current model_number plus this control

      change_initial_net = .true. ! switch nuclear reaction network
      new_net_name = 'pp_and_cno_extras.net'
      
       set_rate_n14pg = 'Imbriani'
       set_rate_c12ag = 'Kunz'

      kappa_file_prefix = 'OP_gs98'
      kappa_lowT_prefix = 'lowT_fa05_gs98' ! for lower temperatures.
     
      change_lnPgas_flag = .true.
      new_lnPgas_flag = .true.

      set_initial_model_number = .true.
      initial_model_number = 0

/ ! end of star_job namelist

&controls

        initial_mass =1.0d0

!asteroseismology
      calculate_Brunt_N2 = .true.
      use_brunt_dlnRho_form = .true.
      num_cells_for_smooth_brunt_B = 5 !2
      use_brunt_gradmuX_form = .true. !.false.
      min_magnitude_brunt_B = 1d-5   !-1d99

      ! check for retries and backups as part of test_suite
      ! you can/should delete this for use outside of test_suite
         max_number_backups = 200
         max_number_retries = 500

!atmosphere boundary conditions
       which_atm_option = 'Eddington_grey'  ! 'Krishna_Swamy' !

      ! parameters for integrate T(tau)
         atm_int_errtol = 1d-7
         dump_int_atm_info_model_number = -1111 ! write atm structure to terminal

      use_Type2_opacities = .false.
      Zbase = 0.02d0

      D_DSI_factor = 0

! iscan for claculating pulsation frequency     
      x_integer_ctrl(5)=800 ! 200
     
!timestep
      max_years_for_timestep = 0.2d8

!mesh adjustment      
       varcontrol_target = 1d-3
      mesh_delta_coeff = 0.3
 
      max_allowed_nz =50000 !8000 ! maximum number of grid points allowed

      min_center_cell_dq = 1d-13
      max_center_cell_dq = 1d-12
      max_surface_cell_dq = 1d-13
      mesh_logX_species(1) = 'h1' ! iso name such as he4
      mesh_logX_min_for_extra(1) = -8
      mesh_min_dr_div_cs = -1 ! limit (in seconds) on sound crossing time for mesh refinement
         ! don't split if sound crossing time would drop below this limit
      xtra_coef_a_l_nb_czb = 0.1 ! above lower nonburn convective boundary
      xtra_coef_b_l_nb_czb = 0.1 ! below lower nonburn convective boundary
      xtra_dist_a_l_nb_czb = 0.01 ! above lower nonburn convective boundary
      xtra_dist_b_l_nb_czb = 0.01 ! below lower nonburn convective boundary
      ! multiply mesh_delta_coeff near convection zone boundary (czb) by the following factors
      ! value < 1 gives increased resolution
     
      ! the center mass fraction of he4 is used to control this extra coefficient
      ! the default settings limit the application to after center he4 is depleted
      xtra_coef_czb_full_on = 0.97 ! if center he4 < this, then use xtra coef's
      xtra_coef_czb_full_off = 0.99 ! if center he4 > this, then don't use xtra coef's
 
      num_cells_for_smooth_brunt_B = 5 !2

      smooth_convective_bdy = .true.
      convective_bdy_weight = 0.01

!when stop
      max_model_number = 2000 ! negative means no maximum
      max_age = 4.57d9 ! in years

      xa_central_lower_limit_species(1) = 'h1'
      xa_central_lower_limit(1) = 0.1

! controls for output
      max_num_profile_models = 2000
      photostep = 100
      profile_interval =100
      history_interval = 1
      terminal_cnt = 10
      write_header_frequency = 100

!stellar winds     
      RGB_wind_scheme = 'Reimers'
      AGB_wind_scheme = 'Blocker'
      RGB_to_AGB_wind_switch = 1d-4
      Reimers_wind_eta = 0.0!0.7d0 
      Blocker_wind_eta = 0.0!0.7d0 

!overshooting     
         ! overshooting depends on the classification of the convective zone
         ! and can be different at the top and the bottom of the zone.
        
      min_overshoot_q = 1d-5
         ! overshooting is only allowed at locations with mass m >= min_overshoot_q * mstar.
         ! e.g., if min_overshoot_q = 0.1, then only the outer 90% by mass can have overshooting.
         ! this provides a simple way of suppressing bogus center overshooting in which a small
         ! convective region at the core can produce excessively large overshooting because of
         ! a large pressure scale height at the center.
           
     D_mix_ov_limit = 1d2
         ! overshooting shuts off when the exponential decay has dropped the Eulerian diffusion
         ! coefficient to this level.
        
         ! parameters for exponential diffusive overshoot are described in the paper by Falk Herwig,
         ! "The evolution of AGB stars with convective overshoot", A&A, 360, 952-968 (2000).
        
         ! NOTE: in addition to giving these 'f' parameters non-zero values, you should also
         ! check the settings for mass_for_overshoot_full_on and mass_for_overshoot_full_off.

     overshoot_f_above_nonburn = 0 !.08
     overshoot_f_below_nonburn = 0 !.08
     overshoot_f_above_burn_h = 0 !.017
     overshoot_f_below_burn_h = 0 !.017
     overshoot_f_above_burn_he = 0
     overshoot_f_below_burn_he = 0
     overshoot_f_above_burn_z = 0
     overshoot_f_below_burn_z = 0
        
     overshoot_below_noburn_factor = 1
         ! multiply overshoot_f_below_nonburn by this factor
         ! during dredge up phase of AGB thermal pulse
     max_DUP_counter = 200
         ! for deciding when to terminate use of overshoot_below_noburn_factor
           
     ovr_below_burn_he_factor = 1
         ! multiply overshoot_f_below_burn_he by this factor
         ! after the first AGB thermal pulse
        
         ! the switch from convective mixing to overshooting happens
         ! at a distance f0*Hp from the estimated location where grad_ad == grad_rad;
         ! where Hp is the pressure scale height at that location.
         ! a value < 0 tells the system to use the value of f for f0.
     overshoot_f0_above_nonburn = -1 !0.2 ! -1
     overshoot_f0_below_nonburn = -1 !0.2 ! -1
     overshoot_f0_above_burn_h = -1 !  -1
     overshoot_f0_below_burn_h = -1 !0.2 !  -1
     overshoot_f0_above_burn_he = -1
     overshoot_f0_below_burn_he = -1
     overshoot_f0_above_burn_z = -1
     overshoot_f0_below_burn_z = -1
                 
        ! optional step function for overshooting to be used in place of exponential
     !overshoot_step_fraction = 0 ! use this if > 0
        
        ! as above, f0*Hp determines r0 where switch from convection to overshooting.
        ! D0 = diffusion coefficient D at point r0
        ! Hp0 = the scale height at r0.
        ! overshooting extends a distance f*Hp0 from r0
        ! with constant diffusion coef D = overshoot_step_fraction*D0
                          
        ! you can specify a range of star masses over which overshooting is gradually enabled
     !mass_for_overshoot_full_on = 1.8  ! 1.8 ! Msun units
     !mass_for_overshoot_full_off = 1.0  ! 1.1 ! Msun units

        ! optional step function for overshooting
        ! note that this can be used simultaneously with exponential overshooting           
     step_overshoot_f_above_nonburn = 0
     step_overshoot_f_below_nonburn = 0
     step_overshoot_f_above_burn_h = 0 !.2
     step_overshoot_f_below_burn_h = 0
     step_overshoot_f_above_burn_he = 0
     step_overshoot_f_below_burn_he = 0
     step_overshoot_f_above_burn_z = 0
     step_overshoot_f_below_burn_z = 0

         ! NOTE: when using step overshoot, you must set f0 as well as f.
                    
     step_overshoot_D = 1d39
     step_overshoot_D0_coeff = 1d60
        
         ! as above, f0*Hp determines r0 where switch from convection to overshooting.
         ! overshooting extends a distance step_f*Hp0 from r0
         ! with constant diffusion coef D = step_D + step_D0_coeff*D0
         ! where D0 = diffusion coefficient D at point r0

!atomic diffusion
      do_element_diffusion =.true. !.false. !   .true. ! determines whether or not we do diffusion
      diffusion_dt_limit = 3.15d10 ! no element diffusion if dt < this limit (in seconds)
      diffusion_T_full_on = 1d3
      diffusion_T_full_off = 1d3
     
      diffusion_calculates_ionization = .true.

      diffusion_use_pure_Coulomb = .true.  !.false.

      ! if false, use atomic diffusion coefficients according to Paquette et al. (1986)


      diffusion_num_classes = 8 ! number of classes of species for diffusion calculations
      diffusion_class_representative(1) = 'h1'
      diffusion_class_representative(2) = 'he3'
      diffusion_class_representative(3) = 'he4'
      diffusion_class_representative(4) = 'c12'
      diffusion_class_representative(5) = 'n14'
      diffusion_class_representative(6) = 'o16'
      diffusion_class_representative(7) = 'ne20'
      diffusion_class_representative(8) = 'mg24'
  
      ! in ascending order.  species goes into 1st class with A_max >= species A
      diffusion_class_A_max(1) = 2
      diffusion_class_A_max(2) = 3
      diffusion_class_A_max(3) = 4
      diffusion_class_A_max(4) = 12
      diffusion_class_A_max(5) = 14
      diffusion_class_A_max(6) = 16
      diffusion_class_A_max(7) = 20
      diffusion_class_A_max(8) = 10000
                
      diffusion_use_isolve = .true.
      diffusion_rtol_for_isolve = 1d-6
      diffusion_atol_for_isolve = 1d-7
      diffusion_maxsteps_for_isolve = 1000
      diffusion_isolve_solver = 'ros2_solver'

      cubic_interpolation_in_Z = .true.
     
/ ! end of controls namelist
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.mesastar.org/pipermail/mesa-users/attachments/20160822/1d4348bd/attachment.html>


More information about the Mesa-users mailing list