SPEED
MAKE_SEISMIC_MOMENT.f90 File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine make_seismic_moment (nn, ct, ww, dd, sxx, syy, szz, syz, szx, sxy, nl_sism, check_ns, check_dist_ns, length_cns, ielem, facsmom, tausmom, func_type, func_indx, func_data, nfdata, nf, t
 Creates the seismic moment tensor.
 

Function/Subroutine Documentation

◆ make_seismic_moment()

subroutine make_seismic_moment ( integer*4  nn,
real*8, dimension(nn)  ct,
real*8, dimension(nn)  ww,
real*8, dimension(nn,nn)  dd,
real*8, dimension(nn,nn,nn)  sxx,
real*8, dimension(nn,nn,nn)  syy,
real*8, dimension(nn,nn,nn)  szz,
real*8, dimension(nn,nn,nn)  syz,
real*8, dimension(nn,nn,nn)  szx,
real*8, dimension(nn,nn,nn)  sxy,
integer*4  nl_sism,
integer*4, dimension(length_cns,4)  check_ns,
real*8, dimension(length_cns,1)  check_dist_ns,
integer*4  length_cns,
integer*4  ielem,
real*8, dimension(nl_sism,6)  facsmom,
real*8, dimension(nl_sism,1)  tausmom,
integer*4, dimension(nf)  func_type,
integer*4, dimension(nf +1)  func_indx,
real*8, dimension(nfdata)  func_data,
integer*4  nfdata,
integer*4  nf,
  t 
)

Creates the seismic moment tensor.

Author
Ilario Mazzieri
Date
September, 2013
Version
1.0
Parameters
[in]nnnumber of 1-D GLL nodes
[in]xqGLL nodes
[in]wqGLL weights
[in]ddspectral derivatives matrix
[in]nl_sismnumber of sismic loads
[in]length_cnslength of check_ns
[in]check_nsinfo about sismic load
[in]check_dist_nsdistance from epicenter
[in]ielemelement index
[in]facsmomseismic moment factor
[in]tausmomrise time factor
[in]numberof fuctions
[in]func_typefunction types
[in]func_indxfunction indixes for function data
[in]func_datafunction data
[in]t_stresstime
[in]cs_nnzlength of cs
[in]csspectral connectivity vector
[in]tag_funcfunctin label
[in]nelem_locnumber of local elements
[in]local_el_numlocal element numeration
[in]nn_locnumber of local nodes
[in]local_n_numlocal node numeration
[in,out]sxxnodal values for the stress tensor
[in,out]syynodal values for the stress tensor
[in,out]szznodal values for the stress tensor
[in,out]syznodal values for the stress tensor
[in,out]szxnodal values for the stress tensor
[in,out]sxynodal values for the stress tensor

Definition at line 53 of file MAKE_SEISMIC_MOMENT.f90.

60 cs_nnz,cs,tag_func, &
61 nelem_loc, local_el_num,&
62 nn_loc, local_n_num)
63
64 implicit none
65
66 integer*4 :: nn, nelem_loc, nn_loc
67 integer*4 :: i,j,k,p,q,r
68 integer*4 :: cs_nnz
69 integer*4 :: is,in,fn,ie
70 integer*4 :: length_cns,ielem,nl_sism
71 integer*4 :: nf,nfdata
72
73 integer*4, dimension(nelem_loc) :: local_el_num
74 integer*4, dimension(nn_loc) :: local_n_num
75 integer*4, dimension(nf) :: func_type
76 integer*4, dimension(nf +1) :: func_indx
77 integer*4, dimension(nf) :: tag_func
78 integer*4, dimension(0:cs_nnz) :: cs
79
80 integer*4, dimension(length_cns,4) :: check_ns
81
82 real*8 :: t_stress
83 real*8 :: get_func_value_sism, val_sism
84
85 real*8, dimension(nn) :: ct,ww
86 real*8, dimension(nfdata) :: func_data
87
88 real*8, dimension(nn,nn) :: dd
89 real*8, dimension(length_cns,1) :: check_dist_ns
90 real*8, dimension(nl_sism,6) :: facsmom
91 real*8, dimension(nl_sism,1) :: tausmom
92
93 real*8, dimension(nn,nn,nn) :: sxx,syy,szz,syz,szx,sxy
94 logical :: exist
95
96 if ((ielem .ge. check_ns(1,4)) .and. (ielem .le. check_ns(length_cns,4))) then
97
98 do i = 1,length_cns
99
100 if (ielem .eq. check_ns(i,4)) then
101
102
103 do r = 1,nn
104 do q = 1,nn
105 do p = 1,nn
106
107 is = nn*nn*(r -1) +nn*(q -1) +p
108
109
110 in = local_n_num(cs(cs(ielem -1) +is))
111
112 do fn = 1,nf
113
114 if ((tag_func(fn) .eq. (check_ns(i,2))) .and. (check_ns(i,1) .eq. in)) then
115
116
117 !STRESS RE-CALCULATION
118
119 !t_stress = current time,
120 !check_dist_ns(i,1) = rupture time
121 !tausmom = rise time
122 !facsmom = seismic moment
123 val_sism = get_func_value_sism(nf,func_type,func_indx,func_data, nfdata, &
124 fn,t_stress,check_dist_ns(i,1), &
125 tausmom(check_ns(i,3),1))
126
127 !write(*,*) VAL_SISM
128 !read(*,*)
129
130! sxx(p,q,r) = sxx(p,q,r) - GET_FUNC_VALUE_SISM(nf,func_type,func_indx,func_data, nfdata, &
131! fn,t_stress,check_dist_ns(i,1), &
132! tausmom(check_ns(i,3),1)) * facsmom(check_ns(i,3),1)
133
134! syy(p,q,r) = syy(p,q,r) - GET_FUNC_VALUE_SISM(nf,func_type,func_indx,func_data, nfdata,&
135! fn,t_stress,check_dist_ns(i,1), &
136! tausmom(check_ns(i,3),1)) * facsmom(check_ns(i,3),2)
137
138! szz(p,q,r) = szz(p,q,r) - GET_FUNC_VALUE_SISM(nf,func_type,func_indx,func_data, nfdata,&
139! fn,t_stress,check_dist_ns(i,1), &
140! tausmom(check_ns(i,3),1)) * facsmom(check_ns(i,3),3)
141
142! syz(p,q,r) = syz(p,q,r) - GET_FUNC_VALUE_SISM(nf,func_type,func_indx,func_data, nfdata,&
143! fn,t_stress,check_dist_ns(i,1), &
144! tausmom(check_ns(i,3),1)) * facsmom(check_ns(i,3),4)
145
146! szx(p,q,r) = szx(p,q,r) - GET_FUNC_VALUE_SISM(nf,func_type,func_indx,func_data, nfdata,&
147! fn,t_stress,check_dist_ns(i,1), &
148! tausmom(check_ns(i,3),1)) * facsmom(check_ns(i,3),5)
149
150! sxy(p,q,r) = sxy(p,q,r) - GET_FUNC_VALUE_SISM(nf,func_type,func_indx,func_data, nfdata,&
151! fn,t_stress,check_dist_ns(i,1), &
152! tausmom(check_ns(i,3),1)) * facsmom(check_ns(i,3),6)
153
154 sxx(p,q,r) = sxx(p,q,r) - val_sism * facsmom(check_ns(i,3),1)
155
156 syy(p,q,r) = syy(p,q,r) - val_sism * facsmom(check_ns(i,3),2)
157
158 szz(p,q,r) = szz(p,q,r) - val_sism * facsmom(check_ns(i,3),3)
159
160 syz(p,q,r) = syz(p,q,r) - val_sism * facsmom(check_ns(i,3),4)
161
162 szx(p,q,r) = szx(p,q,r) - val_sism * facsmom(check_ns(i,3),5)
163
164 sxy(p,q,r) = sxy(p,q,r) - val_sism * facsmom(check_ns(i,3),6)
165
166
167 ! inquire(file="sft_th_debug.txt", exist=exist)
168 ! if (exist) then
169 ! open(82, file="sft_th_debug.txt", status="old", position="append", action="write")
170 ! else
171 ! open(82, file="sft_th_debug.txt", status="new", action="write")
172 ! end if
173 ! write(82,*) 'time', t_stress, 'funct_value', GET_FUNC_VALUE_SISM(nf,func_type,func_indx,&
174 ! func_data, nfdata, fn,t_stress,check_dist_ns(i,1), &
175 ! tausmom(check_ns(i,3),1))
176 ! close(82)
177
178
179 endif ! if ((tag_func(fn).eq.(check_node_sism(i,2))) &...
180
181 enddo !fn = 1,nf
182
183 enddo !p
184 enddo !q
185 enddo !r
186
187 endif !ielem
188
189 enddo !i
190
191
192 endif !ielem
193
194 return
195
real *8 function get_func_value_sism(nb_fnc, type_fnc, ind_fnc, data_fnc, nb_data_fnc, id_fnc,
Computes time evolution function for seismic source.

References get_func_value_sism().

Here is the call graph for this function: