diff --git a/src/constants.f90 b/src/constants.f90 index 69bfefdb04569334a61b26950ef9d6ef7fe24e12..52a98cfbbd549e0af1653294892e3b8631f9180d 100644 --- a/src/constants.f90 +++ b/src/constants.f90 @@ -13,7 +13,7 @@ module constants real(real12), parameter :: hbar = 1.05457148e-34_real12 !Reduced Planck constant real(real12), parameter :: atomic_mass=1.67262158e-27_real12 !Atomic mass of hydrogen real(real12), parameter :: avogadros=6.022e23_real12 !Avogadro's number - real(real12), parameter :: StefBoltz=5.670373e-8_real12 !Stefan-Boltzmann constant + real(real12), parameter :: TINY=1.0e-12_real12 !Tiny number to avoid division by zero integer(int12) , parameter :: fields=6_int12 !Number of fields in the system end module constants diff --git a/src/constructions.f90 b/src/constructions.f90 index ae076aee95b5181952b4308979e728affe8ad291..5550305b96f230bce8f44354664a3ad9bd2a4de8 100644 --- a/src/constructions.f90 +++ b/src/constructions.f90 @@ -27,7 +27,7 @@ module constructions real(real12) :: volume !volume of the block real(real12), dimension(3) :: Length !length of the block in 3 dimensions integer(int12) :: iheater !whether the block is a heater - real(real12) :: kappa, rho, heat_capacity, tau, em !physical properties of the material + real(real12) :: kappa, rho, heat_capacity, tau !physical properties of the material end type heatblock !!Defines the material to be used in the simulation type material @@ -39,7 +39,6 @@ module constructions real(real12) :: rho !density of the material real(real12) :: sound_speed !speed of sound in the material real(real12) :: tau !relaxation time - real(real12) :: em !emissivity of the material logical :: source !whether the material is a source of heat ?? end type material end module constructions diff --git a/src/heating.f90 b/src/heating.f90 index 79ef2c0e9d49c02f27cae39853205ba5a1d2872d..b40d9e2e39c7af6844e04334ebba39bd1e094f8f 100644 --- a/src/heating.f90 +++ b/src/heating.f90 @@ -11,9 +11,9 @@ !!! Author: Harry Mclean, Frank Davies, Steven Hepplestone !!!################################################################################################# module Heating - use constants, only: real12, int12, pi, StefBoltz + use constants, only: real12, int12, pi use globe_data, only: Temp_p, Temp_pp, Heat, heated_volume - use inputs, only: nx,ny,nz, grid, NA, power_in, time_step, T_System, freq, ntime, T_Bath + use inputs, only: nx,ny,nz, grid, NA, power_in, time_step, T_System, freq, ntime use materials, only: material implicit none contains @@ -35,7 +35,7 @@ contains Q = 0._real12 POWER = power_in time = time_step * real(itime,real12) - time_pulse = 1E-12_real12 + time_pulse = 0.5_real12 heated_volume=0.0 heated_num=0 @@ -89,20 +89,20 @@ contains !------------------------------ ! AC oscillatory heating raw, with power correction !------------------------------ - Q(IA) = POWER * (sin(time * 2.0_real12 * PI * freq)**2.0_real12) + Q(IA) = POWER * (sin(time * 2.0_real12 * PI * freq)**2.0_real12)& + +POWER*2.0_real12*PI*freq*tau*sin(2.0_real12*time*2.0_real12*PI*freq) !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ case(5) !------------------------------ - ! AC oscillatory heating raw, with power correction + ! Radiative heating !------------------------------ - Q(IA) = POWER * (sin(time * 2.0_real12 * PI * freq)**2.0_real12)& - +POWER*2.0_real12*PI*freq*tau*sin(2.0_real12*time*2.0_real12*PI*freq) + ! Q(IA)= e*eps* T**4 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ case(6) !------------------------------ - ! Pulsed heating + ! Step one heating !------------------------------ - if ((itime .le. time_pulse)) then + if (itime == 1) then Q(IA) = POWER else Q(IA) = 0.0 @@ -111,15 +111,12 @@ contains end select !------------------------------ - ! If emissitivity is not zero, then calculate the radiative heating + ! convert from Q~density to Q !------------------------------ - - Q(IA) = Q(IA) - grid(ix,iy,iz)%em * grid(ix,iy,iz)%length(1)*& - grid(ix,iy,iz)%length(2)*StefBoltz & - * ((Temp_p(IA)**4.0_real12) - (T_Bath**4.0_real12)) + !Qdens(IA)=Q(IA) + Q(IA)=Q(IA) !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ - !------------------------------ diff --git a/src/inputs.f90 b/src/inputs.f90 index fb40836f8e560973258c1894c4eb1930582f3e1a..160c9968254897f9baf3354a110a8d3dcc361c3f 100644 --- a/src/inputs.f90 +++ b/src/inputs.f90 @@ -525,7 +525,7 @@ subroutine read_mat(unit) type(material), dimension(100) :: dum_mat character(1024) :: buffer integer :: reason, j - integer, dimension(8) :: readvarmat + integer, dimension(7) :: readvarmat integer :: i, index i=0 @@ -571,7 +571,6 @@ subroutine read_mat(unit) CALL assignD(buffer,"rho" ,dum_mat(i)%rho ,readvarmat(5))! assign rho CALL assignD(buffer,"sound_speed" ,dum_mat(i)%sound_speed ,readvarmat(6))! assign sound_speed CALL assignD(buffer,"tau" ,dum_mat(i)%tau ,readvarmat(7))! assign tau - CALL assignD(buffer,"em" ,dum_mat(i)%em ,readvarmat(8))! assign emissivity end do read ! Check for duplicate indices diff --git a/src/material.f90 b/src/material.f90 index 5af32ba2ec974ce0c4f42114480b97afc075e7a6..86e1077c7706aa7c5bc6cf7c0ac06ae23025e186 100644 --- a/src/material.f90 +++ b/src/material.f90 @@ -22,11 +22,11 @@ module materials contains -subroutine material(imaterial_type,kappa,kappa3D,h_conv,heat_capacity,rho,sound_speed,tau, em) +subroutine material(imaterial_type,kappa,kappa3D,h_conv,heat_capacity,rho,sound_speed,tau) integer(int12), intent(in) ::imaterial_type integer(int12) :: i, tmp - real(real12), intent(inout) :: kappa3D, kappa, h_conv, heat_capacity, sound_speed, rho, tau, em + real(real12), intent(inout) :: kappa3D, kappa, h_conv, heat_capacity, sound_speed, rho, tau logical :: found !!!============================================= @@ -66,7 +66,6 @@ subroutine material(imaterial_type,kappa,kappa3D,h_conv,heat_capacity,rho,sound_ rho = input_materials(i)%rho sound_speed = input_materials(i)%sound_speed tau = input_materials(i)%tau - em = input_materials(i)%em exit mat_loop end if diff --git a/src/setup.f90 b/src/setup.f90 index a1591b79887db3c6b327df2e1f0b1e9474b58e22..59d2cd91c26546338b264b238e2f445c792f3bf2 100644 --- a/src/setup.f90 +++ b/src/setup.f90 @@ -32,7 +32,7 @@ module setup !!!################################################################################################# subroutine set_global_variables() integer(int12) :: ix,iy,iz,index - real(real12) :: kappa,kappa3D,h_conv,heat_capacity,rho,sound_speed,tau, em + real(real12) :: kappa,kappa3D,h_conv,heat_capacity,rho,sound_speed,tau allocate(Temp_cur(nx, ny, nz)) allocate(Temp_p(NA)) allocate(Temp_pp(NA)) @@ -51,12 +51,11 @@ module setup do ix = 1, nx index = index + 1 CALL material(grid(ix,iy,iz)%imaterial_type,& - kappa,kappa3D,h_conv,heat_capacity,rho,sound_speed,tau, em) + kappa,kappa3D,h_conv,heat_capacity,rho,sound_speed,tau) grid(ix,iy,iz)%kappa = kappa grid(ix,iy,iz)%rho = rho grid(ix,iy,iz)%heat_capacity = heat_capacity grid(ix,iy,iz)%tau = tau*inverse_time*inverse_time - grid(ix,iy,iz)%em = em lin_rhoc(index) = rho*heat_capacity if (Check_Stability) CALL stability(kappa, rho, heat_capacity, ix, iy, iz) end do