157 ,rg_gll_displacement&
159 ,rg_gll_acceleration&
170 do igll = 1,ig_ngll_total
171 rg_gll_displacement(1,igll) = rg_gll_displacement(1,igll) + rg_dt*rg_gll_velocity(1,igll) + rg_dt2*0.5*rg_gll_acceleration(1,igll)
172 rg_gll_displacement(2,igll) = rg_gll_displacement(2,igll) + rg_dt*rg_gll_velocity(2,igll) + rg_dt2*0.5*rg_gll_acceleration(2,igll)
173 rg_gll_displacement(3,igll) = rg_gll_displacement(3,igll) + rg_dt*rg_gll_velocity(3,igll) + rg_dt2*0.5*rg_gll_acceleration(3,igll)
178 do igll = 1,ig_ngll_total
179 rg_gll_velocity(1,igll) = rg_gll_velocity(1,igll) + rg_dt*(1.0-rg_newmark_gamma)*rg_gll_acceleration(1,igll)
180 rg_gll_velocity(2,igll) = rg_gll_velocity(2,igll) + rg_dt*(1.0-rg_newmark_gamma)*rg_gll_acceleration(2,igll)
181 rg_gll_velocity(3,igll) = rg_gll_velocity(3,igll) + rg_dt*(1.0-rg_newmark_gamma)*rg_gll_acceleration(3,igll)
186 do igll = 1,ig_ngll_total
187 rg_gll_acceleration(1,igll) = 0.0
188 rg_gll_acceleration(2,igll) = 0.0
189 rg_gll_acceleration(3,igll) = 0.0
212 ,rg_gll_acceleration&
220 do igll = 1,ig_ngll_total
221 rg_gll_acceleration(1,igll) = rg_gll_acceleration(1,igll)*rg_gll_mass_matrix(igll)
222 rg_gll_acceleration(2,igll) = rg_gll_acceleration(2,igll)*rg_gll_mass_matrix(igll)
223 rg_gll_acceleration(3,igll) = rg_gll_acceleration(3,igll)*rg_gll_mass_matrix(igll)
226 do igll = 1,ig_ngll_total
227 rg_gll_velocity(1,igll) = rg_gll_velocity(1,igll) + rg_dt*rg_newmark_gamma*rg_gll_acceleration(1,igll)
228 rg_gll_velocity(2,igll) = rg_gll_velocity(2,igll) + rg_dt*rg_newmark_gamma*rg_gll_acceleration(2,igll)
229 rg_gll_velocity(3,igll) = rg_gll_velocity(3,igll) + rg_dt*rg_newmark_gamma*rg_gll_acceleration(3,igll)
253 ,rg_gll_lagrange_deriv&
254 ,rg_gll_displacement&
255 ,rg_gll_acceleration&
260 ,ig_receiver_saving_incr&
272 ,rg_hexa_gll_jacobian_det&
297 integer,
intent(in) :: elt_start
298 integer,
intent(in) :: elt_end
300 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpx1
301 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpx2
302 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpx3
303 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpy1
304 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpy2
305 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpy3
306 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpz1
307 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpz2
308 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpz3
309 real,
dimension(IG_NDOF,IG_NGLL,IG_NGLL,IG_NGLL) :: rl_displacement_gll
310 real,
dimension(IG_NDOF,IG_NGLL,IG_NGLL,IG_NGLL) :: rl_acceleration_gll
374 do iel = elt_start,elt_end
382 igll = ig_hexa_gll_glonum(m,l,k,iel)
384 rl_displacement_gll(1,m,l,k) = rg_gll_displacement(1,igll)
385 rl_displacement_gll(2,m,l,k) = rg_gll_displacement(2,igll)
386 rl_displacement_gll(3,m,l,k) = rg_gll_displacement(3,igll)
403 duxdxi = rl_displacement_gll(1,1,l,k)*rg_gll_lagrange_deriv(1,m) &
404 + rl_displacement_gll(1,2,l,k)*rg_gll_lagrange_deriv(2,m) &
405 + rl_displacement_gll(1,3,l,k)*rg_gll_lagrange_deriv(3,m) &
406 + rl_displacement_gll(1,4,l,k)*rg_gll_lagrange_deriv(4,m) &
407 + rl_displacement_gll(1,5,l,k)*rg_gll_lagrange_deriv(5,m)
409 duydxi = rl_displacement_gll(2,1,l,k)*rg_gll_lagrange_deriv(1,m) &
410 + rl_displacement_gll(2,2,l,k)*rg_gll_lagrange_deriv(2,m) &
411 + rl_displacement_gll(2,3,l,k)*rg_gll_lagrange_deriv(3,m) &
412 + rl_displacement_gll(2,4,l,k)*rg_gll_lagrange_deriv(4,m) &
413 + rl_displacement_gll(2,5,l,k)*rg_gll_lagrange_deriv(5,m)
415 duzdxi = rl_displacement_gll(3,1,l,k)*rg_gll_lagrange_deriv(1,m) &
416 + rl_displacement_gll(3,2,l,k)*rg_gll_lagrange_deriv(2,m) &
417 + rl_displacement_gll(3,3,l,k)*rg_gll_lagrange_deriv(3,m) &
418 + rl_displacement_gll(3,4,l,k)*rg_gll_lagrange_deriv(4,m) &
419 + rl_displacement_gll(3,5,l,k)*rg_gll_lagrange_deriv(5,m)
421 duxdet = rl_displacement_gll(1,m,1,k)*rg_gll_lagrange_deriv(1,l) &
422 + rl_displacement_gll(1,m,2,k)*rg_gll_lagrange_deriv(2,l) &
423 + rl_displacement_gll(1,m,3,k)*rg_gll_lagrange_deriv(3,l) &
424 + rl_displacement_gll(1,m,4,k)*rg_gll_lagrange_deriv(4,l) &
425 + rl_displacement_gll(1,m,5,k)*rg_gll_lagrange_deriv(5,l)
427 duydet = rl_displacement_gll(2,m,1,k)*rg_gll_lagrange_deriv(1,l) &
428 + rl_displacement_gll(2,m,2,k)*rg_gll_lagrange_deriv(2,l) &
429 + rl_displacement_gll(2,m,3,k)*rg_gll_lagrange_deriv(3,l) &
430 + rl_displacement_gll(2,m,4,k)*rg_gll_lagrange_deriv(4,l) &
431 + rl_displacement_gll(2,m,5,k)*rg_gll_lagrange_deriv(5,l)
433 duzdet = rl_displacement_gll(3,m,1,k)*rg_gll_lagrange_deriv(1,l) &
434 + rl_displacement_gll(3,m,2,k)*rg_gll_lagrange_deriv(2,l) &
435 + rl_displacement_gll(3,m,3,k)*rg_gll_lagrange_deriv(3,l) &
436 + rl_displacement_gll(3,m,4,k)*rg_gll_lagrange_deriv(4,l) &
437 + rl_displacement_gll(3,m,5,k)*rg_gll_lagrange_deriv(5,l)
439 duxdze = rl_displacement_gll(1,m,l,1)*rg_gll_lagrange_deriv(1,k) &
440 + rl_displacement_gll(1,m,l,2)*rg_gll_lagrange_deriv(2,k) &
441 + rl_displacement_gll(1,m,l,3)*rg_gll_lagrange_deriv(3,k) &
442 + rl_displacement_gll(1,m,l,4)*rg_gll_lagrange_deriv(4,k) &
443 + rl_displacement_gll(1,m,l,5)*rg_gll_lagrange_deriv(5,k)
445 duydze = rl_displacement_gll(2,m,l,1)*rg_gll_lagrange_deriv(1,k) &
446 + rl_displacement_gll(2,m,l,2)*rg_gll_lagrange_deriv(2,k) &
447 + rl_displacement_gll(2,m,l,3)*rg_gll_lagrange_deriv(3,k) &
448 + rl_displacement_gll(2,m,l,4)*rg_gll_lagrange_deriv(4,k) &
449 + rl_displacement_gll(2,m,l,5)*rg_gll_lagrange_deriv(5,k)
451 duzdze = rl_displacement_gll(3,m,l,1)*rg_gll_lagrange_deriv(1,k) &
452 + rl_displacement_gll(3,m,l,2)*rg_gll_lagrange_deriv(2,k) &
453 + rl_displacement_gll(3,m,l,3)*rg_gll_lagrange_deriv(3,k) &
454 + rl_displacement_gll(3,m,l,4)*rg_gll_lagrange_deriv(4,k) &
455 + rl_displacement_gll(3,m,l,5)*rg_gll_lagrange_deriv(5,k)
459 dxidx = rg_hexa_gll_dxidx(m,l,k,iel)
460 dxidy = rg_hexa_gll_dxidy(m,l,k,iel)
461 dxidz = rg_hexa_gll_dxidz(m,l,k,iel)
462 detdx = rg_hexa_gll_detdx(m,l,k,iel)
463 detdy = rg_hexa_gll_detdy(m,l,k,iel)
464 detdz = rg_hexa_gll_detdz(m,l,k,iel)
465 dzedx = rg_hexa_gll_dzedx(m,l,k,iel)
466 dzedy = rg_hexa_gll_dzedy(m,l,k,iel)
467 dzedz = rg_hexa_gll_dzedz(m,l,k,iel)
469 duxdx = duxdxi*dxidx + duxdet*detdx + duxdze*dzedx
470 duxdy = duxdxi*dxidy + duxdet*detdy + duxdze*dzedy
471 duxdz = duxdxi*dxidz + duxdet*detdz + duxdze*dzedz
472 duydx = duydxi*dxidx + duydet*detdx + duydze*dzedx
473 duydy = duydxi*dxidy + duydet*detdy + duydze*dzedy
474 duydz = duydxi*dxidz + duydet*detdz + duydze*dzedz
475 duzdx = duzdxi*dxidx + duzdet*detdx + duzdze*dzedx
476 duzdy = duzdxi*dxidy + duzdet*detdy + duzdze*dzedy
477 duzdz = duzdxi*dxidz + duzdet*detdz + duzdze*dzedz
480 trace_tau = (rg_hexa_gll_rhovp2(m,l,k,iel) - 2.0*rg_hexa_gll_rhovs2(m,l,k,iel))*(duxdx+duydy+duzdz)
481 tauxx = trace_tau + 2.0*rg_hexa_gll_rhovs2(m,l,k,iel)*duxdx
482 tauyy = trace_tau + 2.0*rg_hexa_gll_rhovs2(m,l,k,iel)*duydy
483 tauzz = trace_tau + 2.0*rg_hexa_gll_rhovs2(m,l,k,iel)*duzdz
484 tauxy = rg_hexa_gll_rhovs2(m,l,k,iel)*(duxdy+duydx)
485 tauxz = rg_hexa_gll_rhovs2(m,l,k,iel)*(duxdz+duzdx)
486 tauyz = rg_hexa_gll_rhovs2(m,l,k,iel)*(duydz+duzdy)
492 do imem_var = 1,ig_nrelax
494 tmpx1 = rg_mem_var_exp(imem_var)
496 tmpx2 = rg_hexa_gll_rhovp2(m,l,k,iel)*rg_hexa_gll_wkqp(imem_var,m,l,k,iel)
497 tmpx3 = 2.0*rg_hexa_gll_rhovs2(m,l,k,iel)*rg_hexa_gll_wkqs(imem_var,m,l,k,iel)
499 tmpx4 = (duxdx+duydy+duzdz)*(tmpx2 - tmpx3)
503 tauxx_n12 = tmpx1*rg_hexa_gll_ksixx(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*duxdx + tmpx4)
504 tauyy_n12 = tmpx1*rg_hexa_gll_ksiyy(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*duydy + tmpx4)
505 tauzz_n12 = tmpx1*rg_hexa_gll_ksizz(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*duzdz + tmpx4)
506 tauxy_n12 = tmpx1*rg_hexa_gll_ksixy(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*0.5*(duxdy+duydx))
507 tauxz_n12 = tmpx1*rg_hexa_gll_ksixz(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*0.5*(duxdz+duzdx))
508 tauyz_n12 = tmpx1*rg_hexa_gll_ksiyz(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*0.5*(duydz+duzdy))
511 tauxx = tauxx - 0.5*(tauxx_n12 + rg_hexa_gll_ksixx(imem_var,m,l,k,iel))
512 tauyy = tauyy - 0.5*(tauyy_n12 + rg_hexa_gll_ksiyy(imem_var,m,l,k,iel))
513 tauzz = tauzz - 0.5*(tauzz_n12 + rg_hexa_gll_ksizz(imem_var,m,l,k,iel))
514 tauxy = tauxy - 0.5*(tauxy_n12 + rg_hexa_gll_ksixy(imem_var,m,l,k,iel))
515 tauxz = tauxz - 0.5*(tauxz_n12 + rg_hexa_gll_ksixz(imem_var,m,l,k,iel))
516 tauyz = tauyz - 0.5*(tauyz_n12 + rg_hexa_gll_ksiyz(imem_var,m,l,k,iel))
519 rg_hexa_gll_ksixx(imem_var,m,l,k,iel) = tauxx_n12
520 rg_hexa_gll_ksiyy(imem_var,m,l,k,iel) = tauyy_n12
521 rg_hexa_gll_ksizz(imem_var,m,l,k,iel) = tauzz_n12
522 rg_hexa_gll_ksixy(imem_var,m,l,k,iel) = tauxy_n12
523 rg_hexa_gll_ksixz(imem_var,m,l,k,iel) = tauxz_n12
524 rg_hexa_gll_ksiyz(imem_var,m,l,k,iel) = tauyz_n12
531 intpx1(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxx*dxidx+tauxy*dxidy+tauxz*dxidz)
532 intpx2(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxx*detdx+tauxy*detdy+tauxz*detdz)
533 intpx3(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxx*dzedx+tauxy*dzedy+tauxz*dzedz)
535 intpy1(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxy*dxidx+tauyy*dxidy+tauyz*dxidz)
536 intpy2(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxy*detdx+tauyy*detdy+tauyz*detdz)
537 intpy3(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxy*dzedx+tauyy*dzedy+tauyz*dzedz)
539 intpz1(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxz*dxidx+tauyz*dxidy+tauzz*dxidz)
540 intpz2(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxz*detdx+tauyz*detdy+tauzz*detdz)
541 intpz3(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxz*dzedx+tauyz*dzedy+tauzz*dzedz)
553 tmpx1 = intpx1(1,l,k)*rg_gll_lagrange_deriv(m,1)*rg_gll_weight(1) &
554 + intpx1(2,l,k)*rg_gll_lagrange_deriv(m,2)*rg_gll_weight(2) &
555 + intpx1(3,l,k)*rg_gll_lagrange_deriv(m,3)*rg_gll_weight(3) &
556 + intpx1(4,l,k)*rg_gll_lagrange_deriv(m,4)*rg_gll_weight(4) &
557 + intpx1(5,l,k)*rg_gll_lagrange_deriv(m,5)*rg_gll_weight(5)
559 tmpy1 = intpy1(1,l,k)*rg_gll_lagrange_deriv(m,1)*rg_gll_weight(1) &
560 + intpy1(2,l,k)*rg_gll_lagrange_deriv(m,2)*rg_gll_weight(2) &
561 + intpy1(3,l,k)*rg_gll_lagrange_deriv(m,3)*rg_gll_weight(3) &
562 + intpy1(4,l,k)*rg_gll_lagrange_deriv(m,4)*rg_gll_weight(4) &
563 + intpy1(5,l,k)*rg_gll_lagrange_deriv(m,5)*rg_gll_weight(5)
565 tmpz1 = intpz1(1,l,k)*rg_gll_lagrange_deriv(m,1)*rg_gll_weight(1) &
566 + intpz1(2,l,k)*rg_gll_lagrange_deriv(m,2)*rg_gll_weight(2) &
567 + intpz1(3,l,k)*rg_gll_lagrange_deriv(m,3)*rg_gll_weight(3) &
568 + intpz1(4,l,k)*rg_gll_lagrange_deriv(m,4)*rg_gll_weight(4) &
569 + intpz1(5,l,k)*rg_gll_lagrange_deriv(m,5)*rg_gll_weight(5)
571 tmpx2 = intpx2(m,1,k)*rg_gll_lagrange_deriv(l,1)*rg_gll_weight(1) &
572 + intpx2(m,2,k)*rg_gll_lagrange_deriv(l,2)*rg_gll_weight(2) &
573 + intpx2(m,3,k)*rg_gll_lagrange_deriv(l,3)*rg_gll_weight(3) &
574 + intpx2(m,4,k)*rg_gll_lagrange_deriv(l,4)*rg_gll_weight(4) &
575 + intpx2(m,5,k)*rg_gll_lagrange_deriv(l,5)*rg_gll_weight(5)
577 tmpy2 = intpy2(m,1,k)*rg_gll_lagrange_deriv(l,1)*rg_gll_weight(1) &
578 + intpy2(m,2,k)*rg_gll_lagrange_deriv(l,2)*rg_gll_weight(2) &
579 + intpy2(m,3,k)*rg_gll_lagrange_deriv(l,3)*rg_gll_weight(3) &
580 + intpy2(m,4,k)*rg_gll_lagrange_deriv(l,4)*rg_gll_weight(4) &
581 + intpy2(m,5,k)*rg_gll_lagrange_deriv(l,5)*rg_gll_weight(5)
583 tmpz2 = intpz2(m,1,k)*rg_gll_lagrange_deriv(l,1)*rg_gll_weight(1) &
584 + intpz2(m,2,k)*rg_gll_lagrange_deriv(l,2)*rg_gll_weight(2) &
585 + intpz2(m,3,k)*rg_gll_lagrange_deriv(l,3)*rg_gll_weight(3) &
586 + intpz2(m,4,k)*rg_gll_lagrange_deriv(l,4)*rg_gll_weight(4) &
587 + intpz2(m,5,k)*rg_gll_lagrange_deriv(l,5)*rg_gll_weight(5)
589 tmpx3 = intpx3(m,l,1)*rg_gll_lagrange_deriv(k,1)*rg_gll_weight(1) &
590 + intpx3(m,l,2)*rg_gll_lagrange_deriv(k,2)*rg_gll_weight(2) &
591 + intpx3(m,l,3)*rg_gll_lagrange_deriv(k,3)*rg_gll_weight(3) &
592 + intpx3(m,l,4)*rg_gll_lagrange_deriv(k,4)*rg_gll_weight(4) &
593 + intpx3(m,l,5)*rg_gll_lagrange_deriv(k,5)*rg_gll_weight(5)
595 tmpy3 = intpy3(m,l,1)*rg_gll_lagrange_deriv(k,1)*rg_gll_weight(1) &
596 + intpy3(m,l,2)*rg_gll_lagrange_deriv(k,2)*rg_gll_weight(2) &
597 + intpy3(m,l,3)*rg_gll_lagrange_deriv(k,3)*rg_gll_weight(3) &
598 + intpy3(m,l,4)*rg_gll_lagrange_deriv(k,4)*rg_gll_weight(4) &
599 + intpy3(m,l,5)*rg_gll_lagrange_deriv(k,5)*rg_gll_weight(5)
601 tmpz3 = intpz3(m,l,1)*rg_gll_lagrange_deriv(k,1)*rg_gll_weight(1) &
602 + intpz3(m,l,2)*rg_gll_lagrange_deriv(k,2)*rg_gll_weight(2) &
603 + intpz3(m,l,3)*rg_gll_lagrange_deriv(k,3)*rg_gll_weight(3) &
604 + intpz3(m,l,4)*rg_gll_lagrange_deriv(k,4)*rg_gll_weight(4) &
605 + intpz3(m,l,5)*rg_gll_lagrange_deriv(k,5)*rg_gll_weight(5)
607 fac1 = rg_gll_weight(l)*rg_gll_weight(k)
608 fac2 = rg_gll_weight(m)*rg_gll_weight(k)
609 fac3 = rg_gll_weight(m)*rg_gll_weight(l)
611 rl_acceleration_gll(1,m,l,k) = (fac1*tmpx1 + fac2*tmpx2 + fac3*tmpx3)
612 rl_acceleration_gll(2,m,l,k) = (fac1*tmpy1 + fac2*tmpy2 + fac3*tmpy3)
613 rl_acceleration_gll(3,m,l,k) = (fac1*tmpz1 + fac2*tmpz2 + fac3*tmpz3)
623 igll = ig_hexa_gll_glonum(m,l,k,iel)
625 rg_gll_acceleration(1,igll) = rg_gll_acceleration(1,igll) - rl_acceleration_gll(1,m,l,k)
626 rg_gll_acceleration(2,igll) = rg_gll_acceleration(2,igll) - rl_acceleration_gll(2,m,l,k)
627 rg_gll_acceleration(3,igll) = rg_gll_acceleration(3,igll) - rl_acceleration_gll(3,m,l,k)
656 ,rg_gll_lagrange_deriv&
657 ,rg_gll_displacement&
658 ,rg_gll_acceleration&
663 ,ig_receiver_saving_incr&
675 ,rg_hexa_gll_jacobian_det&
700 integer,
intent(in) :: elt_start
701 integer,
intent(in) :: elt_end
703 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpx1
704 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpx2
705 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpx3
706 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpy1
707 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpy2
708 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpy3
709 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpz1
710 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpz2
711 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpz3
712 real,
dimension(IG_NDOF,IG_NGLL,IG_NGLL,IG_NGLL) :: rl_displacement_gll
713 real,
dimension(IG_NDOF,IG_NGLL,IG_NGLL,IG_NGLL) :: rl_acceleration_gll
777 do iel = elt_start,elt_end
785 rl_acceleration_gll(1,m,l,k) = 0.0
786 rl_acceleration_gll(2,m,l,k) = 0.0
787 rl_acceleration_gll(3,m,l,k) = 0.0
799 igll = ig_hexa_gll_glonum(m,l,k,iel)
801 rl_displacement_gll(1,m,l,k) = rg_gll_displacement(1,igll)
802 rl_displacement_gll(2,m,l,k) = rg_gll_displacement(2,igll)
803 rl_displacement_gll(3,m,l,k) = rg_gll_displacement(3,igll)
820 duxdxi = rl_displacement_gll(1,1,l,k)*rg_gll_lagrange_deriv(1,m) &
821 + rl_displacement_gll(1,2,l,k)*rg_gll_lagrange_deriv(2,m) &
822 + rl_displacement_gll(1,3,l,k)*rg_gll_lagrange_deriv(3,m) &
823 + rl_displacement_gll(1,4,l,k)*rg_gll_lagrange_deriv(4,m) &
824 + rl_displacement_gll(1,5,l,k)*rg_gll_lagrange_deriv(5,m) &
825 + rl_displacement_gll(1,6,l,k)*rg_gll_lagrange_deriv(6,m)
827 duydxi = rl_displacement_gll(2,1,l,k)*rg_gll_lagrange_deriv(1,m) &
828 + rl_displacement_gll(2,2,l,k)*rg_gll_lagrange_deriv(2,m) &
829 + rl_displacement_gll(2,3,l,k)*rg_gll_lagrange_deriv(3,m) &
830 + rl_displacement_gll(2,4,l,k)*rg_gll_lagrange_deriv(4,m) &
831 + rl_displacement_gll(2,5,l,k)*rg_gll_lagrange_deriv(5,m) &
832 + rl_displacement_gll(2,6,l,k)*rg_gll_lagrange_deriv(6,m)
834 duzdxi = rl_displacement_gll(3,1,l,k)*rg_gll_lagrange_deriv(1,m) &
835 + rl_displacement_gll(3,2,l,k)*rg_gll_lagrange_deriv(2,m) &
836 + rl_displacement_gll(3,3,l,k)*rg_gll_lagrange_deriv(3,m) &
837 + rl_displacement_gll(3,4,l,k)*rg_gll_lagrange_deriv(4,m) &
838 + rl_displacement_gll(3,5,l,k)*rg_gll_lagrange_deriv(5,m) &
839 + rl_displacement_gll(3,6,l,k)*rg_gll_lagrange_deriv(6,m)
841 duxdet = rl_displacement_gll(1,m,1,k)*rg_gll_lagrange_deriv(1,l) &
842 + rl_displacement_gll(1,m,2,k)*rg_gll_lagrange_deriv(2,l) &
843 + rl_displacement_gll(1,m,3,k)*rg_gll_lagrange_deriv(3,l) &
844 + rl_displacement_gll(1,m,4,k)*rg_gll_lagrange_deriv(4,l) &
845 + rl_displacement_gll(1,m,5,k)*rg_gll_lagrange_deriv(5,l) &
846 + rl_displacement_gll(1,m,6,k)*rg_gll_lagrange_deriv(6,l)
848 duydet = rl_displacement_gll(2,m,1,k)*rg_gll_lagrange_deriv(1,l) &
849 + rl_displacement_gll(2,m,2,k)*rg_gll_lagrange_deriv(2,l) &
850 + rl_displacement_gll(2,m,3,k)*rg_gll_lagrange_deriv(3,l) &
851 + rl_displacement_gll(2,m,4,k)*rg_gll_lagrange_deriv(4,l) &
852 + rl_displacement_gll(2,m,5,k)*rg_gll_lagrange_deriv(5,l) &
853 + rl_displacement_gll(2,m,6,k)*rg_gll_lagrange_deriv(6,l)
855 duzdet = rl_displacement_gll(3,m,1,k)*rg_gll_lagrange_deriv(1,l) &
856 + rl_displacement_gll(3,m,2,k)*rg_gll_lagrange_deriv(2,l) &
857 + rl_displacement_gll(3,m,3,k)*rg_gll_lagrange_deriv(3,l) &
858 + rl_displacement_gll(3,m,4,k)*rg_gll_lagrange_deriv(4,l) &
859 + rl_displacement_gll(3,m,5,k)*rg_gll_lagrange_deriv(5,l) &
860 + rl_displacement_gll(3,m,6,k)*rg_gll_lagrange_deriv(6,l)
862 duxdze = rl_displacement_gll(1,m,l,1)*rg_gll_lagrange_deriv(1,k) &
863 + rl_displacement_gll(1,m,l,2)*rg_gll_lagrange_deriv(2,k) &
864 + rl_displacement_gll(1,m,l,3)*rg_gll_lagrange_deriv(3,k) &
865 + rl_displacement_gll(1,m,l,4)*rg_gll_lagrange_deriv(4,k) &
866 + rl_displacement_gll(1,m,l,5)*rg_gll_lagrange_deriv(5,k) &
867 + rl_displacement_gll(1,m,l,6)*rg_gll_lagrange_deriv(6,k)
869 duydze = rl_displacement_gll(2,m,l,1)*rg_gll_lagrange_deriv(1,k) &
870 + rl_displacement_gll(2,m,l,2)*rg_gll_lagrange_deriv(2,k) &
871 + rl_displacement_gll(2,m,l,3)*rg_gll_lagrange_deriv(3,k) &
872 + rl_displacement_gll(2,m,l,4)*rg_gll_lagrange_deriv(4,k) &
873 + rl_displacement_gll(2,m,l,5)*rg_gll_lagrange_deriv(5,k) &
874 + rl_displacement_gll(2,m,l,6)*rg_gll_lagrange_deriv(6,k)
876 duzdze = rl_displacement_gll(3,m,l,1)*rg_gll_lagrange_deriv(1,k) &
877 + rl_displacement_gll(3,m,l,2)*rg_gll_lagrange_deriv(2,k) &
878 + rl_displacement_gll(3,m,l,3)*rg_gll_lagrange_deriv(3,k) &
879 + rl_displacement_gll(3,m,l,4)*rg_gll_lagrange_deriv(4,k) &
880 + rl_displacement_gll(3,m,l,5)*rg_gll_lagrange_deriv(5,k) &
881 + rl_displacement_gll(3,m,l,6)*rg_gll_lagrange_deriv(6,k)
885 dxidx = rg_hexa_gll_dxidx(m,l,k,iel)
886 dxidy = rg_hexa_gll_dxidy(m,l,k,iel)
887 dxidz = rg_hexa_gll_dxidz(m,l,k,iel)
888 detdx = rg_hexa_gll_detdx(m,l,k,iel)
889 detdy = rg_hexa_gll_detdy(m,l,k,iel)
890 detdz = rg_hexa_gll_detdz(m,l,k,iel)
891 dzedx = rg_hexa_gll_dzedx(m,l,k,iel)
892 dzedy = rg_hexa_gll_dzedy(m,l,k,iel)
893 dzedz = rg_hexa_gll_dzedz(m,l,k,iel)
895 duxdx = duxdxi*dxidx + duxdet*detdx + duxdze*dzedx
896 duxdy = duxdxi*dxidy + duxdet*detdy + duxdze*dzedy
897 duxdz = duxdxi*dxidz + duxdet*detdz + duxdze*dzedz
898 duydx = duydxi*dxidx + duydet*detdx + duydze*dzedx
899 duydy = duydxi*dxidy + duydet*detdy + duydze*dzedy
900 duydz = duydxi*dxidz + duydet*detdz + duydze*dzedz
901 duzdx = duzdxi*dxidx + duzdet*detdx + duzdze*dzedx
902 duzdy = duzdxi*dxidy + duzdet*detdy + duzdze*dzedy
903 duzdz = duzdxi*dxidz + duzdet*detdz + duzdze*dzedz
906 trace_tau = (rg_hexa_gll_rhovp2(m,l,k,iel) - 2.0*rg_hexa_gll_rhovs2(m,l,k,iel))*(duxdx+duydy+duzdz)
907 tauxx = trace_tau + 2.0*rg_hexa_gll_rhovs2(m,l,k,iel)*duxdx
908 tauyy = trace_tau + 2.0*rg_hexa_gll_rhovs2(m,l,k,iel)*duydy
909 tauzz = trace_tau + 2.0*rg_hexa_gll_rhovs2(m,l,k,iel)*duzdz
910 tauxy = rg_hexa_gll_rhovs2(m,l,k,iel)*(duxdy+duydx)
911 tauxz = rg_hexa_gll_rhovs2(m,l,k,iel)*(duxdz+duzdx)
912 tauyz = rg_hexa_gll_rhovs2(m,l,k,iel)*(duydz+duzdy)
918 do imem_var = 1,ig_nrelax
920 tmpx1 = rg_mem_var_exp(imem_var)
922 tmpx2 = rg_hexa_gll_rhovp2(m,l,k,iel)*rg_hexa_gll_wkqp(imem_var,m,l,k,iel)
923 tmpx3 = 2.0*rg_hexa_gll_rhovs2(m,l,k,iel)*rg_hexa_gll_wkqs(imem_var,m,l,k,iel)
925 tmpx4 = (duxdx+duydy+duzdz)*(tmpx2 - tmpx3)
929 tauxx_n12 = tmpx1*rg_hexa_gll_ksixx(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*duxdx + tmpx4)
930 tauyy_n12 = tmpx1*rg_hexa_gll_ksiyy(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*duydy + tmpx4)
931 tauzz_n12 = tmpx1*rg_hexa_gll_ksizz(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*duzdz + tmpx4)
932 tauxy_n12 = tmpx1*rg_hexa_gll_ksixy(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*0.5*(duxdy+duydx))
933 tauxz_n12 = tmpx1*rg_hexa_gll_ksixz(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*0.5*(duxdz+duzdx))
934 tauyz_n12 = tmpx1*rg_hexa_gll_ksiyz(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*0.5*(duydz+duzdy))
937 tauxx = tauxx - 0.5*(tauxx_n12 + rg_hexa_gll_ksixx(imem_var,m,l,k,iel))
938 tauyy = tauyy - 0.5*(tauyy_n12 + rg_hexa_gll_ksiyy(imem_var,m,l,k,iel))
939 tauzz = tauzz - 0.5*(tauzz_n12 + rg_hexa_gll_ksizz(imem_var,m,l,k,iel))
940 tauxy = tauxy - 0.5*(tauxy_n12 + rg_hexa_gll_ksixy(imem_var,m,l,k,iel))
941 tauxz = tauxz - 0.5*(tauxz_n12 + rg_hexa_gll_ksixz(imem_var,m,l,k,iel))
942 tauyz = tauyz - 0.5*(tauyz_n12 + rg_hexa_gll_ksiyz(imem_var,m,l,k,iel))
945 rg_hexa_gll_ksixx(imem_var,m,l,k,iel) = tauxx_n12
946 rg_hexa_gll_ksiyy(imem_var,m,l,k,iel) = tauyy_n12
947 rg_hexa_gll_ksizz(imem_var,m,l,k,iel) = tauzz_n12
948 rg_hexa_gll_ksixy(imem_var,m,l,k,iel) = tauxy_n12
949 rg_hexa_gll_ksixz(imem_var,m,l,k,iel) = tauxz_n12
950 rg_hexa_gll_ksiyz(imem_var,m,l,k,iel) = tauyz_n12
957 intpx1(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxx*dxidx+tauxy*dxidy+tauxz*dxidz)
958 intpx2(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxx*detdx+tauxy*detdy+tauxz*detdz)
959 intpx3(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxx*dzedx+tauxy*dzedy+tauxz*dzedz)
961 intpy1(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxy*dxidx+tauyy*dxidy+tauyz*dxidz)
962 intpy2(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxy*detdx+tauyy*detdy+tauyz*detdz)
963 intpy3(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxy*dzedx+tauyy*dzedy+tauyz*dzedz)
965 intpz1(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxz*dxidx+tauyz*dxidy+tauzz*dxidz)
966 intpz2(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxz*detdx+tauyz*detdy+tauzz*detdz)
967 intpz3(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxz*dzedx+tauyz*dzedy+tauzz*dzedz)
978 tmpx1 = intpx1(1,l,k)*rg_gll_lagrange_deriv(m,1)*rg_gll_weight(1) &
979 + intpx1(2,l,k)*rg_gll_lagrange_deriv(m,2)*rg_gll_weight(2) &
980 + intpx1(3,l,k)*rg_gll_lagrange_deriv(m,3)*rg_gll_weight(3) &
981 + intpx1(4,l,k)*rg_gll_lagrange_deriv(m,4)*rg_gll_weight(4) &
982 + intpx1(5,l,k)*rg_gll_lagrange_deriv(m,5)*rg_gll_weight(5) &
983 + intpx1(6,l,k)*rg_gll_lagrange_deriv(m,6)*rg_gll_weight(6)
985 tmpy1 = intpy1(1,l,k)*rg_gll_lagrange_deriv(m,1)*rg_gll_weight(1) &
986 + intpy1(2,l,k)*rg_gll_lagrange_deriv(m,2)*rg_gll_weight(2) &
987 + intpy1(3,l,k)*rg_gll_lagrange_deriv(m,3)*rg_gll_weight(3) &
988 + intpy1(4,l,k)*rg_gll_lagrange_deriv(m,4)*rg_gll_weight(4) &
989 + intpy1(5,l,k)*rg_gll_lagrange_deriv(m,5)*rg_gll_weight(5) &
990 + intpy1(6,l,k)*rg_gll_lagrange_deriv(m,6)*rg_gll_weight(6)
992 tmpz1 = intpz1(1,l,k)*rg_gll_lagrange_deriv(m,1)*rg_gll_weight(1) &
993 + intpz1(2,l,k)*rg_gll_lagrange_deriv(m,2)*rg_gll_weight(2) &
994 + intpz1(3,l,k)*rg_gll_lagrange_deriv(m,3)*rg_gll_weight(3) &
995 + intpz1(4,l,k)*rg_gll_lagrange_deriv(m,4)*rg_gll_weight(4) &
996 + intpz1(5,l,k)*rg_gll_lagrange_deriv(m,5)*rg_gll_weight(5) &
997 + intpz1(6,l,k)*rg_gll_lagrange_deriv(m,6)*rg_gll_weight(6)
999 tmpx2 = intpx2(m,1,k)*rg_gll_lagrange_deriv(l,1)*rg_gll_weight(1) &
1000 + intpx2(m,2,k)*rg_gll_lagrange_deriv(l,2)*rg_gll_weight(2) &
1001 + intpx2(m,3,k)*rg_gll_lagrange_deriv(l,3)*rg_gll_weight(3) &
1002 + intpx2(m,4,k)*rg_gll_lagrange_deriv(l,4)*rg_gll_weight(4) &
1003 + intpx2(m,5,k)*rg_gll_lagrange_deriv(l,5)*rg_gll_weight(5) &
1004 + intpx2(m,6,k)*rg_gll_lagrange_deriv(l,6)*rg_gll_weight(6)
1006 tmpy2 = intpy2(m,1,k)*rg_gll_lagrange_deriv(l,1)*rg_gll_weight(1) &
1007 + intpy2(m,2,k)*rg_gll_lagrange_deriv(l,2)*rg_gll_weight(2) &
1008 + intpy2(m,3,k)*rg_gll_lagrange_deriv(l,3)*rg_gll_weight(3) &
1009 + intpy2(m,4,k)*rg_gll_lagrange_deriv(l,4)*rg_gll_weight(4) &
1010 + intpy2(m,5,k)*rg_gll_lagrange_deriv(l,5)*rg_gll_weight(5) &
1011 + intpy2(m,6,k)*rg_gll_lagrange_deriv(l,6)*rg_gll_weight(6)
1013 tmpz2 = intpz2(m,1,k)*rg_gll_lagrange_deriv(l,1)*rg_gll_weight(1) &
1014 + intpz2(m,2,k)*rg_gll_lagrange_deriv(l,2)*rg_gll_weight(2) &
1015 + intpz2(m,3,k)*rg_gll_lagrange_deriv(l,3)*rg_gll_weight(3) &
1016 + intpz2(m,4,k)*rg_gll_lagrange_deriv(l,4)*rg_gll_weight(4) &
1017 + intpz2(m,5,k)*rg_gll_lagrange_deriv(l,5)*rg_gll_weight(5) &
1018 + intpz2(m,6,k)*rg_gll_lagrange_deriv(l,6)*rg_gll_weight(6)
1020 tmpx3 = intpx3(m,l,1)*rg_gll_lagrange_deriv(k,1)*rg_gll_weight(1) &
1021 + intpx3(m,l,2)*rg_gll_lagrange_deriv(k,2)*rg_gll_weight(2) &
1022 + intpx3(m,l,3)*rg_gll_lagrange_deriv(k,3)*rg_gll_weight(3) &
1023 + intpx3(m,l,4)*rg_gll_lagrange_deriv(k,4)*rg_gll_weight(4) &
1024 + intpx3(m,l,5)*rg_gll_lagrange_deriv(k,5)*rg_gll_weight(5) &
1025 + intpx3(m,l,6)*rg_gll_lagrange_deriv(k,6)*rg_gll_weight(6)
1027 tmpy3 = intpy3(m,l,1)*rg_gll_lagrange_deriv(k,1)*rg_gll_weight(1) &
1028 + intpy3(m,l,2)*rg_gll_lagrange_deriv(k,2)*rg_gll_weight(2) &
1029 + intpy3(m,l,3)*rg_gll_lagrange_deriv(k,3)*rg_gll_weight(3) &
1030 + intpy3(m,l,4)*rg_gll_lagrange_deriv(k,4)*rg_gll_weight(4) &
1031 + intpy3(m,l,5)*rg_gll_lagrange_deriv(k,5)*rg_gll_weight(5) &
1032 + intpy3(m,l,6)*rg_gll_lagrange_deriv(k,6)*rg_gll_weight(6)
1034 tmpz3 = intpz3(m,l,1)*rg_gll_lagrange_deriv(k,1)*rg_gll_weight(1) &
1035 + intpz3(m,l,2)*rg_gll_lagrange_deriv(k,2)*rg_gll_weight(2) &
1036 + intpz3(m,l,3)*rg_gll_lagrange_deriv(k,3)*rg_gll_weight(3) &
1037 + intpz3(m,l,4)*rg_gll_lagrange_deriv(k,4)*rg_gll_weight(4) &
1038 + intpz3(m,l,5)*rg_gll_lagrange_deriv(k,5)*rg_gll_weight(5) &
1039 + intpz3(m,l,6)*rg_gll_lagrange_deriv(k,6)*rg_gll_weight(6)
1041 fac1 = rg_gll_weight(l)*rg_gll_weight(k)
1042 fac2 = rg_gll_weight(m)*rg_gll_weight(k)
1043 fac3 = rg_gll_weight(m)*rg_gll_weight(l)
1045 rl_acceleration_gll(1,m,l,k) = rl_acceleration_gll(1,m,l,k) + (fac1*tmpx1 + fac2*tmpx2 + fac3*tmpx3)
1046 rl_acceleration_gll(2,m,l,k) = rl_acceleration_gll(2,m,l,k) + (fac1*tmpy1 + fac2*tmpy2 + fac3*tmpy3)
1047 rl_acceleration_gll(3,m,l,k) = rl_acceleration_gll(3,m,l,k) + (fac1*tmpz1 + fac2*tmpz2 + fac3*tmpz3)
1057 igll = ig_hexa_gll_glonum(m,l,k,iel)
1059 rg_gll_acceleration(1,igll) = rg_gll_acceleration(1,igll) - rl_acceleration_gll(1,m,l,k)
1060 rg_gll_acceleration(2,igll) = rg_gll_acceleration(2,igll) - rl_acceleration_gll(2,m,l,k)
1061 rg_gll_acceleration(3,igll) = rg_gll_acceleration(3,igll) - rl_acceleration_gll(3,m,l,k)
1090 ,rg_gll_lagrange_deriv&
1091 ,rg_gll_displacement&
1092 ,rg_gll_acceleration&
1097 ,ig_receiver_saving_incr&
1099 ,ig_hexa_gll_glonum&
1109 ,rg_hexa_gll_jacobian_det&
1113 ,rg_hexa_gll_rhovs2&
1114 ,rg_hexa_gll_rhovp2&
1115 ,rg_hexa_gll_rhovs2&
1116 ,rg_hexa_gll_rhovp2&
1134 integer,
intent(in) :: elt_start
1135 integer,
intent(in) :: elt_end
1137 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpx1
1138 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpx2
1139 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpx3
1140 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpy1
1141 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpy2
1142 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpy3
1143 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpz1
1144 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpz2
1145 real,
dimension(IG_NGLL,IG_NGLL,IG_NGLL) :: intpz3
1146 real,
dimension(IG_NDOF,IG_NGLL,IG_NGLL,IG_NGLL) :: rl_displacement_gll
1147 real,
dimension(IG_NDOF,IG_NGLL,IG_NGLL,IG_NGLL) :: rl_acceleration_gll
1210 do iel = elt_start,elt_end
1218 rl_acceleration_gll(1,m,l,k) = 0.0
1219 rl_acceleration_gll(2,m,l,k) = 0.0
1220 rl_acceleration_gll(3,m,l,k) = 0.0
1232 igll = ig_hexa_gll_glonum(m,l,k,iel)
1234 rl_displacement_gll(1,m,l,k) = rg_gll_displacement(1,igll)
1235 rl_displacement_gll(2,m,l,k) = rg_gll_displacement(2,igll)
1236 rl_displacement_gll(3,m,l,k) = rg_gll_displacement(3,igll)
1253 duxdxi = rl_displacement_gll(1,1,l,k)*rg_gll_lagrange_deriv(1,m) &
1254 + rl_displacement_gll(1,2,l,k)*rg_gll_lagrange_deriv(2,m) &
1255 + rl_displacement_gll(1,3,l,k)*rg_gll_lagrange_deriv(3,m) &
1256 + rl_displacement_gll(1,4,l,k)*rg_gll_lagrange_deriv(4,m) &
1257 + rl_displacement_gll(1,5,l,k)*rg_gll_lagrange_deriv(5,m) &
1258 + rl_displacement_gll(1,6,l,k)*rg_gll_lagrange_deriv(6,m) &
1259 + rl_displacement_gll(1,7,l,k)*rg_gll_lagrange_deriv(7,m)
1261 duydxi = rl_displacement_gll(2,1,l,k)*rg_gll_lagrange_deriv(1,m) &
1262 + rl_displacement_gll(2,2,l,k)*rg_gll_lagrange_deriv(2,m) &
1263 + rl_displacement_gll(2,3,l,k)*rg_gll_lagrange_deriv(3,m) &
1264 + rl_displacement_gll(2,4,l,k)*rg_gll_lagrange_deriv(4,m) &
1265 + rl_displacement_gll(2,5,l,k)*rg_gll_lagrange_deriv(5,m) &
1266 + rl_displacement_gll(2,6,l,k)*rg_gll_lagrange_deriv(6,m) &
1267 + rl_displacement_gll(2,7,l,k)*rg_gll_lagrange_deriv(7,m)
1269 duzdxi = rl_displacement_gll(3,1,l,k)*rg_gll_lagrange_deriv(1,m) &
1270 + rl_displacement_gll(3,2,l,k)*rg_gll_lagrange_deriv(2,m) &
1271 + rl_displacement_gll(3,3,l,k)*rg_gll_lagrange_deriv(3,m) &
1272 + rl_displacement_gll(3,4,l,k)*rg_gll_lagrange_deriv(4,m) &
1273 + rl_displacement_gll(3,5,l,k)*rg_gll_lagrange_deriv(5,m) &
1274 + rl_displacement_gll(3,6,l,k)*rg_gll_lagrange_deriv(6,m) &
1275 + rl_displacement_gll(3,7,l,k)*rg_gll_lagrange_deriv(7,m)
1277 duxdet = rl_displacement_gll(1,m,1,k)*rg_gll_lagrange_deriv(1,l) &
1278 + rl_displacement_gll(1,m,2,k)*rg_gll_lagrange_deriv(2,l) &
1279 + rl_displacement_gll(1,m,3,k)*rg_gll_lagrange_deriv(3,l) &
1280 + rl_displacement_gll(1,m,4,k)*rg_gll_lagrange_deriv(4,l) &
1281 + rl_displacement_gll(1,m,5,k)*rg_gll_lagrange_deriv(5,l) &
1282 + rl_displacement_gll(1,m,6,k)*rg_gll_lagrange_deriv(6,l) &
1283 + rl_displacement_gll(1,m,7,k)*rg_gll_lagrange_deriv(7,l)
1285 duydet = rl_displacement_gll(2,m,1,k)*rg_gll_lagrange_deriv(1,l) &
1286 + rl_displacement_gll(2,m,2,k)*rg_gll_lagrange_deriv(2,l) &
1287 + rl_displacement_gll(2,m,3,k)*rg_gll_lagrange_deriv(3,l) &
1288 + rl_displacement_gll(2,m,4,k)*rg_gll_lagrange_deriv(4,l) &
1289 + rl_displacement_gll(2,m,5,k)*rg_gll_lagrange_deriv(5,l) &
1290 + rl_displacement_gll(2,m,6,k)*rg_gll_lagrange_deriv(6,l) &
1291 + rl_displacement_gll(2,m,7,k)*rg_gll_lagrange_deriv(7,l)
1293 duzdet = rl_displacement_gll(3,m,1,k)*rg_gll_lagrange_deriv(1,l) &
1294 + rl_displacement_gll(3,m,2,k)*rg_gll_lagrange_deriv(2,l) &
1295 + rl_displacement_gll(3,m,3,k)*rg_gll_lagrange_deriv(3,l) &
1296 + rl_displacement_gll(3,m,4,k)*rg_gll_lagrange_deriv(4,l) &
1297 + rl_displacement_gll(3,m,5,k)*rg_gll_lagrange_deriv(5,l) &
1298 + rl_displacement_gll(3,m,6,k)*rg_gll_lagrange_deriv(6,l) &
1299 + rl_displacement_gll(3,m,7,k)*rg_gll_lagrange_deriv(7,l)
1301 duxdze = rl_displacement_gll(1,m,l,1)*rg_gll_lagrange_deriv(1,k) &
1302 + rl_displacement_gll(1,m,l,2)*rg_gll_lagrange_deriv(2,k) &
1303 + rl_displacement_gll(1,m,l,3)*rg_gll_lagrange_deriv(3,k) &
1304 + rl_displacement_gll(1,m,l,4)*rg_gll_lagrange_deriv(4,k) &
1305 + rl_displacement_gll(1,m,l,5)*rg_gll_lagrange_deriv(5,k) &
1306 + rl_displacement_gll(1,m,l,6)*rg_gll_lagrange_deriv(6,k) &
1307 + rl_displacement_gll(1,m,l,7)*rg_gll_lagrange_deriv(7,k)
1309 duydze = rl_displacement_gll(2,m,l,1)*rg_gll_lagrange_deriv(1,k) &
1310 + rl_displacement_gll(2,m,l,2)*rg_gll_lagrange_deriv(2,k) &
1311 + rl_displacement_gll(2,m,l,3)*rg_gll_lagrange_deriv(3,k) &
1312 + rl_displacement_gll(2,m,l,4)*rg_gll_lagrange_deriv(4,k) &
1313 + rl_displacement_gll(2,m,l,5)*rg_gll_lagrange_deriv(5,k) &
1314 + rl_displacement_gll(2,m,l,6)*rg_gll_lagrange_deriv(6,k) &
1315 + rl_displacement_gll(2,m,l,7)*rg_gll_lagrange_deriv(7,k)
1317 duzdze = rl_displacement_gll(3,m,l,1)*rg_gll_lagrange_deriv(1,k) &
1318 + rl_displacement_gll(3,m,l,2)*rg_gll_lagrange_deriv(2,k) &
1319 + rl_displacement_gll(3,m,l,3)*rg_gll_lagrange_deriv(3,k) &
1320 + rl_displacement_gll(3,m,l,4)*rg_gll_lagrange_deriv(4,k) &
1321 + rl_displacement_gll(3,m,l,5)*rg_gll_lagrange_deriv(5,k) &
1322 + rl_displacement_gll(3,m,l,6)*rg_gll_lagrange_deriv(6,k) &
1323 + rl_displacement_gll(3,m,l,7)*rg_gll_lagrange_deriv(7,k)
1327 dxidx = rg_hexa_gll_dxidx(m,l,k,iel)
1328 dxidy = rg_hexa_gll_dxidy(m,l,k,iel)
1329 dxidz = rg_hexa_gll_dxidz(m,l,k,iel)
1330 detdx = rg_hexa_gll_detdx(m,l,k,iel)
1331 detdy = rg_hexa_gll_detdy(m,l,k,iel)
1332 detdz = rg_hexa_gll_detdz(m,l,k,iel)
1333 dzedx = rg_hexa_gll_dzedx(m,l,k,iel)
1334 dzedy = rg_hexa_gll_dzedy(m,l,k,iel)
1335 dzedz = rg_hexa_gll_dzedz(m,l,k,iel)
1337 duxdx = duxdxi*dxidx + duxdet*detdx + duxdze*dzedx
1338 duxdy = duxdxi*dxidy + duxdet*detdy + duxdze*dzedy
1339 duxdz = duxdxi*dxidz + duxdet*detdz + duxdze*dzedz
1340 duydx = duydxi*dxidx + duydet*detdx + duydze*dzedx
1341 duydy = duydxi*dxidy + duydet*detdy + duydze*dzedy
1342 duydz = duydxi*dxidz + duydet*detdz + duydze*dzedz
1343 duzdx = duzdxi*dxidx + duzdet*detdx + duzdze*dzedx
1344 duzdy = duzdxi*dxidy + duzdet*detdy + duzdze*dzedy
1345 duzdz = duzdxi*dxidz + duzdet*detdz + duzdze*dzedz
1348 trace_tau = (rg_hexa_gll_rhovp2(m,l,k,iel) - 2.0*rg_hexa_gll_rhovs2(m,l,k,iel))*(duxdx+duydy+duzdz)
1349 tauxx = trace_tau + 2.0*rg_hexa_gll_rhovs2(m,l,k,iel)*duxdx
1350 tauyy = trace_tau + 2.0*rg_hexa_gll_rhovs2(m,l,k,iel)*duydy
1351 tauzz = trace_tau + 2.0*rg_hexa_gll_rhovs2(m,l,k,iel)*duzdz
1352 tauxy = rg_hexa_gll_rhovs2(m,l,k,iel)*(duxdy+duydx)
1353 tauxz = rg_hexa_gll_rhovs2(m,l,k,iel)*(duxdz+duzdx)
1354 tauyz = rg_hexa_gll_rhovs2(m,l,k,iel)*(duydz+duzdy)
1360 do imem_var = 1,ig_nrelax
1362 tmpx1 = rg_mem_var_exp(imem_var)
1364 tmpx2 = rg_hexa_gll_rhovp2(m,l,k,iel)*rg_hexa_gll_wkqp(imem_var,m,l,k,iel)
1365 tmpx3 = 2.0*rg_hexa_gll_rhovs2(m,l,k,iel)*rg_hexa_gll_wkqs(imem_var,m,l,k,iel)
1367 tmpx4 = (duxdx+duydy+duzdz)*(tmpx2 - tmpx3)
1371 tauxx_n12 = tmpx1*rg_hexa_gll_ksixx(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*duxdx + tmpx4)
1372 tauyy_n12 = tmpx1*rg_hexa_gll_ksiyy(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*duydy + tmpx4)
1373 tauzz_n12 = tmpx1*rg_hexa_gll_ksizz(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*duzdz + tmpx4)
1374 tauxy_n12 = tmpx1*rg_hexa_gll_ksixy(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*0.5*(duxdy+duydx))
1375 tauxz_n12 = tmpx1*rg_hexa_gll_ksixz(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*0.5*(duxdz+duzdx))
1376 tauyz_n12 = tmpx1*rg_hexa_gll_ksiyz(imem_var,m,l,k,iel) + (1.0 - tmpx1) * (tmpx3*0.5*(duydz+duzdy))
1379 tauxx = tauxx - 0.5*(tauxx_n12 + rg_hexa_gll_ksixx(imem_var,m,l,k,iel))
1380 tauyy = tauyy - 0.5*(tauyy_n12 + rg_hexa_gll_ksiyy(imem_var,m,l,k,iel))
1381 tauzz = tauzz - 0.5*(tauzz_n12 + rg_hexa_gll_ksizz(imem_var,m,l,k,iel))
1382 tauxy = tauxy - 0.5*(tauxy_n12 + rg_hexa_gll_ksixy(imem_var,m,l,k,iel))
1383 tauxz = tauxz - 0.5*(tauxz_n12 + rg_hexa_gll_ksixz(imem_var,m,l,k,iel))
1384 tauyz = tauyz - 0.5*(tauyz_n12 + rg_hexa_gll_ksiyz(imem_var,m,l,k,iel))
1387 rg_hexa_gll_ksixx(imem_var,m,l,k,iel) = tauxx_n12
1388 rg_hexa_gll_ksiyy(imem_var,m,l,k,iel) = tauyy_n12
1389 rg_hexa_gll_ksizz(imem_var,m,l,k,iel) = tauzz_n12
1390 rg_hexa_gll_ksixy(imem_var,m,l,k,iel) = tauxy_n12
1391 rg_hexa_gll_ksixz(imem_var,m,l,k,iel) = tauxz_n12
1392 rg_hexa_gll_ksiyz(imem_var,m,l,k,iel) = tauyz_n12
1399 intpx1(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxx*dxidx+tauxy*dxidy+tauxz*dxidz)
1400 intpx2(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxx*detdx+tauxy*detdy+tauxz*detdz)
1401 intpx3(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxx*dzedx+tauxy*dzedy+tauxz*dzedz)
1403 intpy1(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxy*dxidx+tauyy*dxidy+tauyz*dxidz)
1404 intpy2(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxy*detdx+tauyy*detdy+tauyz*detdz)
1405 intpy3(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxy*dzedx+tauyy*dzedy+tauyz*dzedz)
1407 intpz1(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxz*dxidx+tauyz*dxidy+tauzz*dxidz)
1408 intpz2(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxz*detdx+tauyz*detdy+tauzz*detdz)
1409 intpz3(m,l,k) = rg_hexa_gll_jacobian_det(m,l,k,iel)*(tauxz*dzedx+tauyz*dzedy+tauzz*dzedz)
1420 tmpx1 = intpx1(1,l,k)*rg_gll_lagrange_deriv(m,1)*rg_gll_weight(1) &
1421 + intpx1(2,l,k)*rg_gll_lagrange_deriv(m,2)*rg_gll_weight(2) &
1422 + intpx1(3,l,k)*rg_gll_lagrange_deriv(m,3)*rg_gll_weight(3) &
1423 + intpx1(4,l,k)*rg_gll_lagrange_deriv(m,4)*rg_gll_weight(4) &
1424 + intpx1(5,l,k)*rg_gll_lagrange_deriv(m,5)*rg_gll_weight(5) &
1425 + intpx1(6,l,k)*rg_gll_lagrange_deriv(m,6)*rg_gll_weight(6) &
1426 + intpx1(7,l,k)*rg_gll_lagrange_deriv(m,7)*rg_gll_weight(7)
1428 tmpy1 = intpy1(1,l,k)*rg_gll_lagrange_deriv(m,1)*rg_gll_weight(1) &
1429 + intpy1(2,l,k)*rg_gll_lagrange_deriv(m,2)*rg_gll_weight(2) &
1430 + intpy1(3,l,k)*rg_gll_lagrange_deriv(m,3)*rg_gll_weight(3) &
1431 + intpy1(4,l,k)*rg_gll_lagrange_deriv(m,4)*rg_gll_weight(4) &
1432 + intpy1(5,l,k)*rg_gll_lagrange_deriv(m,5)*rg_gll_weight(5) &
1433 + intpy1(6,l,k)*rg_gll_lagrange_deriv(m,6)*rg_gll_weight(6) &
1434 + intpy1(7,l,k)*rg_gll_lagrange_deriv(m,7)*rg_gll_weight(7)
1436 tmpz1 = intpz1(1,l,k)*rg_gll_lagrange_deriv(m,1)*rg_gll_weight(1) &
1437 + intpz1(2,l,k)*rg_gll_lagrange_deriv(m,2)*rg_gll_weight(2) &
1438 + intpz1(3,l,k)*rg_gll_lagrange_deriv(m,3)*rg_gll_weight(3) &
1439 + intpz1(4,l,k)*rg_gll_lagrange_deriv(m,4)*rg_gll_weight(4) &
1440 + intpz1(5,l,k)*rg_gll_lagrange_deriv(m,5)*rg_gll_weight(5) &
1441 + intpz1(6,l,k)*rg_gll_lagrange_deriv(m,6)*rg_gll_weight(6) &
1442 + intpz1(7,l,k)*rg_gll_lagrange_deriv(m,7)*rg_gll_weight(7)
1444 tmpx2 = intpx2(m,1,k)*rg_gll_lagrange_deriv(l,1)*rg_gll_weight(1) &
1445 + intpx2(m,2,k)*rg_gll_lagrange_deriv(l,2)*rg_gll_weight(2) &
1446 + intpx2(m,3,k)*rg_gll_lagrange_deriv(l,3)*rg_gll_weight(3) &
1447 + intpx2(m,4,k)*rg_gll_lagrange_deriv(l,4)*rg_gll_weight(4) &
1448 + intpx2(m,5,k)*rg_gll_lagrange_deriv(l,5)*rg_gll_weight(5) &
1449 + intpx2(m,6,k)*rg_gll_lagrange_deriv(l,6)*rg_gll_weight(6) &
1450 + intpx2(m,7,k)*rg_gll_lagrange_deriv(l,7)*rg_gll_weight(7)
1452 tmpy2 = intpy2(m,1,k)*rg_gll_lagrange_deriv(l,1)*rg_gll_weight(1) &
1453 + intpy2(m,2,k)*rg_gll_lagrange_deriv(l,2)*rg_gll_weight(2) &
1454 + intpy2(m,3,k)*rg_gll_lagrange_deriv(l,3)*rg_gll_weight(3) &
1455 + intpy2(m,4,k)*rg_gll_lagrange_deriv(l,4)*rg_gll_weight(4) &
1456 + intpy2(m,5,k)*rg_gll_lagrange_deriv(l,5)*rg_gll_weight(5) &
1457 + intpy2(m,6,k)*rg_gll_lagrange_deriv(l,6)*rg_gll_weight(6) &
1458 + intpy2(m,7,k)*rg_gll_lagrange_deriv(l,7)*rg_gll_weight(7)
1460 tmpz2 = intpz2(m,1,k)*rg_gll_lagrange_deriv(l,1)*rg_gll_weight(1) &
1461 + intpz2(m,2,k)*rg_gll_lagrange_deriv(l,2)*rg_gll_weight(2) &
1462 + intpz2(m,3,k)*rg_gll_lagrange_deriv(l,3)*rg_gll_weight(3) &
1463 + intpz2(m,4,k)*rg_gll_lagrange_deriv(l,4)*rg_gll_weight(4) &
1464 + intpz2(m,5,k)*rg_gll_lagrange_deriv(l,5)*rg_gll_weight(5) &
1465 + intpz2(m,6,k)*rg_gll_lagrange_deriv(l,6)*rg_gll_weight(6) &
1466 + intpz2(m,7,k)*rg_gll_lagrange_deriv(l,7)*rg_gll_weight(7)
1468 tmpx3 = intpx3(m,l,1)*rg_gll_lagrange_deriv(k,1)*rg_gll_weight(1) &
1469 + intpx3(m,l,2)*rg_gll_lagrange_deriv(k,2)*rg_gll_weight(2) &
1470 + intpx3(m,l,3)*rg_gll_lagrange_deriv(k,3)*rg_gll_weight(3) &
1471 + intpx3(m,l,4)*rg_gll_lagrange_deriv(k,4)*rg_gll_weight(4) &
1472 + intpx3(m,l,5)*rg_gll_lagrange_deriv(k,5)*rg_gll_weight(5) &
1473 + intpx3(m,l,6)*rg_gll_lagrange_deriv(k,6)*rg_gll_weight(6) &
1474 + intpx3(m,l,7)*rg_gll_lagrange_deriv(k,7)*rg_gll_weight(7)
1476 tmpy3 = intpy3(m,l,1)*rg_gll_lagrange_deriv(k,1)*rg_gll_weight(1) &
1477 + intpy3(m,l,2)*rg_gll_lagrange_deriv(k,2)*rg_gll_weight(2) &
1478 + intpy3(m,l,3)*rg_gll_lagrange_deriv(k,3)*rg_gll_weight(3) &
1479 + intpy3(m,l,4)*rg_gll_lagrange_deriv(k,4)*rg_gll_weight(4) &
1480 + intpy3(m,l,5)*rg_gll_lagrange_deriv(k,5)*rg_gll_weight(5) &
1481 + intpy3(m,l,6)*rg_gll_lagrange_deriv(k,6)*rg_gll_weight(6) &
1482 + intpy3(m,l,7)*rg_gll_lagrange_deriv(k,7)*rg_gll_weight(7)
1484 tmpz3 = intpz3(m,l,1)*rg_gll_lagrange_deriv(k,1)*rg_gll_weight(1) &
1485 + intpz3(m,l,2)*rg_gll_lagrange_deriv(k,2)*rg_gll_weight(2) &
1486 + intpz3(m,l,3)*rg_gll_lagrange_deriv(k,3)*rg_gll_weight(3) &
1487 + intpz3(m,l,4)*rg_gll_lagrange_deriv(k,4)*rg_gll_weight(4) &
1488 + intpz3(m,l,5)*rg_gll_lagrange_deriv(k,5)*rg_gll_weight(5) &
1489 + intpz3(m,l,6)*rg_gll_lagrange_deriv(k,6)*rg_gll_weight(6) &
1490 + intpz3(m,l,7)*rg_gll_lagrange_deriv(k,7)*rg_gll_weight(7)
1492 fac1 = rg_gll_weight(l)*rg_gll_weight(k)
1493 fac2 = rg_gll_weight(m)*rg_gll_weight(k)
1494 fac3 = rg_gll_weight(m)*rg_gll_weight(l)
1496 rl_acceleration_gll(1,m,l,k) = rl_acceleration_gll(1,m,l,k) + (fac1*tmpx1 + fac2*tmpx2 + fac3*tmpx3)
1497 rl_acceleration_gll(2,m,l,k) = rl_acceleration_gll(2,m,l,k) + (fac1*tmpy1 + fac2*tmpy2 + fac3*tmpy3)
1498 rl_acceleration_gll(3,m,l,k) = rl_acceleration_gll(3,m,l,k) + (fac1*tmpz1 + fac2*tmpz2 + fac3*tmpz3)
1508 igll = ig_hexa_gll_glonum(m,l,k,iel)
1510 rg_gll_acceleration(1,igll) = rg_gll_acceleration(1,igll) - rl_acceleration_gll(1,m,l,k)
1511 rg_gll_acceleration(2,igll) = rg_gll_acceleration(2,igll) - rl_acceleration_gll(2,m,l,k)
1512 rg_gll_acceleration(3,igll) = rg_gll_acceleration(3,igll) - rl_acceleration_gll(3,m,l,k)
1538 ,rg_gll_acceleration&
1544 ,rg_quadp_gll_jaco_det&
1545 ,rg_quadp_gll_normal&
1546 ,ig_quadp_gll_glonum&
1547 ,rg_quadp_gll_rhovp&
1569 do iquad = 1,ig_nquad_parax
1574 jaco = rg_quadp_gll_jaco_det(l,k,iquad)
1575 nx = rg_quadp_gll_normal(1,l,k,iquad)
1576 ny = rg_quadp_gll_normal(2,l,k,iquad)
1577 nz = rg_quadp_gll_normal(3,l,k,iquad)
1579 igll = ig_quadp_gll_glonum(l,k,iquad)
1581 vx = rg_gll_velocity(1,igll)
1582 vy = rg_gll_velocity(2,igll)
1583 vz = rg_gll_velocity(3,igll)
1585 vn = vx*nx+vy*ny+vz*nz
1587 tx = rg_quadp_gll_rhovp(l,k,iquad)*vn*nx + rg_quadp_gll_rhovs(l,k,iquad)*(vx-vn*nx)
1588 ty = rg_quadp_gll_rhovp(l,k,iquad)*vn*ny + rg_quadp_gll_rhovs(l,k,iquad)*(vy-vn*ny)
1589 tz = rg_quadp_gll_rhovp(l,k,iquad)*vn*nz + rg_quadp_gll_rhovs(l,k,iquad)*(vz-vn*nz)
1591 rg_gll_acceleration(1,igll) = rg_gll_acceleration(1,igll) - jaco*rg_gll_weight(l)*rg_gll_weight(k)*tx
1592 rg_gll_acceleration(2,igll) = rg_gll_acceleration(2,igll) - jaco*rg_gll_weight(l)*rg_gll_weight(k)*ty
1593 rg_gll_acceleration(3,igll) = rg_gll_acceleration(3,igll) - jaco*rg_gll_weight(l)*rg_gll_weight(k)*tz
1614 ,rg_simu_current_time&
1615 ,rg_gll_acceleration&
1620 ,rg_dcsource_gll_force&
1621 ,ig_hexa_gll_glonum&
1624 ,rg_dcsource_user_func&
1625 ,rg_sfsource_user_func
1634 real,
save :: val_dc
1635 real,
save :: val_pf
1651 do iso = 1,ig_ndcsource
1653 if (tg_dcsource(iso)%icur == 3)
then
1655 val =
gabor(rg_simu_current_time,tg_dcsource(iso)%shift_time,tg_dcsource(iso)%rise_time,1.0)
1657 elseif (tg_dcsource(iso)%icur == 4)
then
1659 val =
expcos(rg_simu_current_time,tg_dcsource(iso)%shift_time,tg_dcsource(iso)%rise_time)
1661 elseif (tg_dcsource(iso)%icur == 5)
then
1663 if (ig_idt == 1) val_dc = 0.0
1664 call
ispli3(rg_simu_current_time,tg_dcsource(iso)%rise_time,s)
1665 val_dc = val_dc + s*rg_dt
1668 elseif (tg_dcsource(iso)%icur == 6)
then
1670 val =
ricker(rg_simu_current_time,tg_dcsource(iso)%shift_time,tg_dcsource(iso)%rise_time,1.0)
1672 elseif (tg_dcsource(iso)%icur == 7)
then
1674 val =
spiexp(rg_simu_current_time,tg_dcsource(iso)%rise_time,1.0)
1676 elseif (tg_dcsource(iso)%icur == 8)
then
1678 val =
fctanh(rg_simu_current_time,tg_dcsource(iso)%shift_time,tg_dcsource(iso)%rise_time,1.0)
1680 elseif (tg_dcsource(iso)%icur == 9)
then
1682 val =
fctanh_dt(rg_simu_current_time,tg_dcsource(iso)%shift_time,tg_dcsource(iso)%rise_time,1.0)
1684 elseif (tg_dcsource(iso)%icur == 10)
then
1686 val = rg_dcsource_user_func(ig_idt)
1693 igll = ig_hexa_gll_glonum(m,l,k,tg_dcsource(iso)%iel)
1694 rg_gll_acceleration(1,igll) = rg_gll_acceleration(1,igll) + rg_dcsource_gll_force(1,m,l,k,iso)*val
1695 rg_gll_acceleration(2,igll) = rg_gll_acceleration(2,igll) + rg_dcsource_gll_force(2,m,l,k,iso)*val
1696 rg_gll_acceleration(3,igll) = rg_gll_acceleration(3,igll) + rg_dcsource_gll_force(3,m,l,k,iso)*val
1709 do ipf = 1,ig_nsfsource
1711 if (tg_sfsource(ipf)%icur == 3)
then
1713 val =
gabor(rg_simu_current_time,tg_sfsource(ipf)%shift_time,tg_sfsource(ipf)%rise_time,1.0)
1715 elseif (tg_sfsource(ipf)%icur == 4)
then
1717 val =
expcos(rg_simu_current_time,tg_sfsource(ipf)%shift_time,tg_sfsource(ipf)%rise_time)
1719 elseif (tg_sfsource(ipf)%icur == 5)
then
1721 if (ig_idt == 1) val_pf = 0.0
1722 call
ispli3(rg_simu_current_time,tg_sfsource(ipf)%rise_time,s)
1723 val_pf = val_pf + s*rg_dt
1726 elseif (tg_sfsource(ipf)%icur == 6)
then
1728 val =
ricker(rg_simu_current_time,tg_sfsource(ipf)%shift_time,tg_sfsource(ipf)%rise_time,1.0)
1730 elseif (tg_sfsource(ipf)%icur == 7)
then
1732 val =
spiexp(rg_simu_current_time,tg_sfsource(ipf)%rise_time,1.0)
1734 elseif (tg_sfsource(ipf)%icur == 8)
then
1736 val =
fctanh(rg_simu_current_time,tg_sfsource(ipf)%shift_time,tg_sfsource(ipf)%rise_time,1.0)
1738 elseif (tg_sfsource(ipf)%icur == 9)
then
1740 val =
fctanh_dt(rg_simu_current_time,tg_sfsource(ipf)%shift_time,tg_sfsource(ipf)%rise_time,1.0)
1742 elseif (tg_sfsource(ipf)%icur == 10)
then
1744 val = rg_sfsource_user_func(ig_idt)
1748 igll = tg_sfsource(ipf)%iequ
1749 idir = tg_sfsource(ipf)%idir
1750 fac = tg_sfsource(ipf)%fac
1751 rg_gll_acceleration(idir,igll) = rg_gll_acceleration(idir,igll) + val*fac
subroutine, public compute_external_force()
This subroutine sets external forces of the system for double couple and single force point sources...
subroutine, public compute_internal_forces_order5(elt_start, elt_end)
This subroutine computes internal forces for spectral-elements of order 5. Stress-strain relationshi...
This module contains subroutines to compute Newmark explicit time marching scheme, external forces , internal forces and boundary traction forces of the system .
real function, public spiexp(t, tp, a)
function to compute spice project exponential source function :
real function, public gabor(t, ts, l, a)
function to compute real part of Gabor wavelet :
This module defines all global variables of EFISPEC3D. Scalar variables are initialized directly in t...
subroutine, public compute_internal_forces_order4(elt_start, elt_end)
This subroutine computes internal forces for spectral-elements of order 4. Stress-strain relationshi...
real function, public ricker(t, ts, tp, a)
function to compute order 2 Ricker wavelet :
real function, public expcos(t, ts, fp)
function to compute euroseistest project source function for case can1 :
This module contains functions and subroutines to compute tsource functions's time history...
subroutine, public compute_absorption_forces()
This subroutine computes absorption forces for any spectral-elements order. A so-called 'P1' explici...
subroutine, public ispli3(t, du, s)
function to compute order 3 spline : see Wikipedia
real function, public fctanh_dt(t, ts, tp, a)
function to compute first derivative of hyperbolic tangent function : with
real function, public fctanh(t, ts, tp, a)
function to compute hyperbolic tangent function (also called 'Bouchon pulse') : ...
subroutine, public compute_internal_forces_order6(elt_start, elt_end)
This subroutine computes internal forces for spectral-elements of order 6. Stress-strain relationshi...
subroutine, public newmark_end()
This subroutine finalizes Newmark time marching scheme at step n+1.
subroutine, public newmark_ini()
This subroutine initializes Newmark time marching scheme at step n+1.