[mesa-users] solar model calibration
Earl Bellinger
bellinger at mps.mpg.de
Mon Aug 22 11:09:54 EDT 2016
Dear Tao,
Basu & Antia 1997 use the method of helioseismic inversions to measure
R_cz from observations of the Sun. The evolutionary models they use as
reference for performing the inversions do not match the observed R_cz
either. Evolutionary models are not yet able to fully model the
acoustic structure of the Sun, including R_cz, so it is expected that
your model does not reproduce R_cz. If having an evolutionary model
with the correct R_cz is important to your study, you can try to add
that constraint to your optimization procedure to penalize deviations
from the measured R_cz, just like you have done with age, R, L, and
Z/X. However, you still may not land on the right values. You could
perhaps also modify the structure of your model to make it match the
structure of the Sun -- this then ceases to be an evolutionary model,
but it depends on what you want.
I hope this answer helps.
Best regards,
Earl Bellinger
Stellar Ages & Galactic Evolution Group
Max-Planck-Institut für Sonnensystemforschung
On Mon, Aug 22, 2016 at 7:38 AM, WU T. <ggwutao at 126.com> wrote:
>
> 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
>
>
>
>
>
> ------------------------------------------------------------------------------
>
> _______________________________________________
> mesa-users mailing list
> mesa-users at lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/mesa-users
>
More information about the Mesa-users
mailing list