SPEED
MAKE_EXPL_SOURCE.f90
Go to the documentation of this file.
1! Copyright (C) 2012 The SPEED FOUNDATION
2! Author: Ilario Mazzieri
3!
4! This file is part of SPEED.
5!
6! SPEED is free software; you can redistribute it and/or modify it
7! under the terms of the GNU Affero General Public License as
8! published by the Free Software Foundation, either version 3 of the
9! License, or (at your option) any later version.
10!
11! SPEED is distributed in the hope that it will be useful, but
12! WITHOUT ANY WARRANTY; without even the implied warranty of
13! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14! Affero General Public License for more details.
15!
16! You should have received a copy of the GNU Affero General Public License
17! along with SPEED. If not, see <http://www.gnu.org/licenses/>.
18
51
52
53 subroutine make_expl_source(nn,ct,ww,dd,&
54 sxx,syy,szz,syz,szx,sxy,&
55 nl_expl,&
56 check_ne,check_dist_ne,&
57 length_cne,ielem,facsexpl,&
58 func_type,func_indx,func_data,nfdata,nf,t_stress,&
59 cs_nnz,cs,tag_func, &
60 nelem_loc, local_el_num,&
61 nn_loc, local_n_num)
62
63 implicit none
64
65 integer*4 :: nn, nelem_loc, nn_loc
66 integer*4 :: i,j,k,p,q,r
67 integer*4 :: length_cne,ielem,nl_expl,nfdata
68 integer*4 :: nf
69 integer*4 :: is,in,fn, ie
70 integer*4 :: cs_nnz
71
72 integer*4, dimension(nelem_loc) :: local_el_num
73 integer*4, dimension(nn_loc) :: local_n_num
74 integer*4, dimension(nf) :: tag_func
75 integer*4, dimension(nf) :: func_type
76 integer*4, dimension(nf+1) :: func_indx
77
78 integer*4, dimension(0:cs_nnz) :: cs
79 integer*4, dimension(length_cne,4) :: check_ne
80
81 real*8 :: get_func_value_sism
82 real*8 :: t_stress
83
84 real*8, dimension(nfdata) :: func_data
85 real*8, dimension(nn) :: ct,ww
86
87 real*8, dimension(nn,nn) :: dd
88 real*8, dimension(length_cne,1) :: check_dist_ne
89 real*8, dimension(nl_expl,6) :: facsexpl
90
91 real*8, dimension(nn,nn,nn) :: sxx,syy,szz,syz,szx,sxy
92
93! STRESS CALCULATION
94
95
96
97 if ((ielem.ge.check_ne(1,4)).and. (ielem.le.check_ne(length_cne,4))) then
98
99 do fn = 1,nf
100
101 do i = 1,length_cne
102 if (ielem.eq.check_ne(i,4)) then
103
104
105 do r = 1,nn
106 do q = 1,nn
107 do p = 1,nn
108
109
110 in = local_n_num(cs(cs(ielem -1) +is))
111
112 if ((tag_func(fn).eq.(check_ne(i,2))) .and.(check_ne(i,1).eq.in)) then
113
114
115 !STRESS RE-CALCULATION
116
117
118 sxx(p,q,r) = sxx(p,q,r) &
119 - get_func_value_sism(nf,func_type,func_indx,func_data, nfdata,&
120 fn,t_stress,check_dist_ne(i,1),0.d0) &
121 * facsexpl(check_ne(i,3),1)
122
123 syy(p,q,r) = syy(p,q,r) &
124 - get_func_value_sism(nf,func_type,func_indx,func_data, nfdata,&
125 fn,t_stress,check_dist_ne(i,1),0.d0) &
126 * facsexpl(check_ne(i,3),2)
127
128 szz(p,q,r) = szz(p,q,r) &
129 - get_func_value_sism(nf,func_type,func_indx,func_data, nfdata,&
130 fn,t_stress,check_dist_ne(i,1),0.d0) &
131 * facsexpl(check_ne(i,3),3)
132
133 syz(p,q,r) = syz(p,q,r) &
134 - get_func_value_sism(nf,func_type,func_indx,func_data, nfdata,&
135 fn,t_stress,check_dist_ne(i,1),0.d0) &
136 * facsexpl(check_ne(i,3),4)
137
138 szx(p,q,r) = szx(p,q,r) &
139 - get_func_value_sism(nf,func_type,func_indx,func_data, nfdata,&
140 fn,t_stress,check_dist_ne(i,1),0.d0) &
141 * facsexpl(check_ne(i,3),5)
142
143 sxy(p,q,r) = sxy(p,q,r) &
144 - get_func_value_sism(nf,func_type,func_indx,func_data, nfdata,&
145 fn,t_stress,check_dist_ne(i,1),0.d0) &
146 * facsexpl(check_ne(i,3),6)
147
148 endif ! if ((tag_func(fn).eq.(check_node_expl(i,2))) &...
149
150 enddo
151 enddo
152 enddo
153
154 endif !ielem
155 enddo !i
156 enddo !fn
157
158 endif !ielem
159
160
161 return
162
163 end subroutine make_expl_source
164
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.
subroutine make_expl_source(nn, ct, ww, dd, sxx, syy, szz, syz, szx, sxy, nl_expl, check_ne, check_dist_ne, length_cne, ielem, facsexpl, func_type, func_indx, func_data, nfdata, nf, t
Computes the explosive moment tensor.