****************************** 
* 
* Carl W. Gable
* gable@lanl.gov
* LaGriT Control File
* http://lagrit.lanl.gov
* http://meshing.lanl.gov
* 
****************************** 
* Add an integer node attribute vector (if_fault)
* Add an real    node attribute vector (slip)
*
cmo / addatt / CMO_NAME / if_fault / vint    / scalar / nnodes
cmo / addatt / CMO_NAME / slip     / vdouble / scalar / nnodes
*
* Make a pset of the fault nodes.
*
define / ZBOT / -16001.0
pset / p_fault / geom/xyz/ 1 0 0 /11999.9 0. ZBOT/12000.01 16000.01 1.0
*
* Set nodes of attribute if_fault to:
* 0 = not fault node
* 1 = fault node
*
cmo / setatt / CMO_NAME / if_fault / 1 0 0 / 0
cmo / setatt / CMO_NAME / slip / 1 0 0 / 0.0
cmo / setatt / CMO_NAME / if_fault / pset get p_fault / 1
cmo / setatt / CMO_NAME / slip     / pset get p_fault / 1.0
*
* Slip needs to taper from 1.0 at 12km to 0.0 at 16km
*
* Find the tapering nodes
*
pset / py12 / attribute / yic / pset get p_fault /  12000. / gt
pset / pz12 / attribute / zic / pset get p_fault / -12000. / lt 
pset / ptaper / union / py12 pz12
pset / p_yz   / inter / py12 pz12
pset / p_y    / not   / py12 pz12
pset / p_z    / not   / pz12 py12
*
* 
cmo / setatt / CMO_NAME / if_fault / pset get ptaper / 2
cmo / setatt / CMO_NAME / if_fault / pset get p_yz / 3
cmo / setatt / CMO_NAME / if_fault / pset get p_y / 4
cmo / setatt / CMO_NAME / if_fault / pset get p_z / 5
*
* 
* Nodes in the p_y list have slip = a*y + b
* a = -0.00025, b = 4
* y = 12000  slip = -3 + 4 = 1
* y = 16000  slip = -4 + 4 = 0
*
* Nodes in the p_z list have slip = a*z + b
* a = 0.00025, b = 4
* z = -12000  slip = -3 + 4 = 1
* y = -16000  slip = -4 + 4 = 0
*
* Nodes in the p_yz list
* Apply both criteria above for y and z and then take the smaller of the two
*
math / multiply / CMO_NAME / slip / pset get p_y / CMO_NAME / yic / -0.00025
math / multiply / CMO_NAME / slip / pset get p_z / CMO_NAME / zic /  0.00025
math / add      / CMO_NAME / slip / pset get p_y / CMO_NAME / slip / 4.0
math / add      / CMO_NAME / slip / pset get p_z / CMO_NAME / slip / 4.0

cmo / addatt / CMO_NAME / tmpvar1 / vdouble / scalar / nnodes
cmo / addatt / CMO_NAME / tmpvar2 / vdouble / scalar / nnodes
math / multiply / CMO_NAME / tmpvar1 / pset get p_yz / CMO_NAME / yic / -0.00025
math / multiply / CMO_NAME / tmpvar2 / pset get p_yz / CMO_NAME / zic /  0.00025
math / add      / CMO_NAME / tmpvar1 / pset get p_yz / CMO_NAME / tmpvar1 / 4.0
math / add      / CMO_NAME / tmpvar2 / pset get p_yz / CMO_NAME / tmpvar2 / 4.0
math / min      / CMO_NAME / slip    / pset get p_yz / &
                  CMO_NAME / tmpvar1 / CMO_NAME / tmpvar2 /
cmo / DELATT / CMO_NAME / tmpvar1
cmo / DELATT / CMO_NAME / tmpvar2
*
* Due to numerical roundoff, a few nodes get set to -1.e-12. Set anything negative to zero.
*
pset / pzero / attribute / slip / 1 0 0 /  0.0 / lt
cmo / setatt / CMO_NAME / slip / pset get pzero / 0.0
*
finish