| 
            Last change
 on this file since 17989 was             17989, checked in by Mathieu Morlighem, 11 years ago           | 
        
        
          | 
             
merged trunk-jpl and trunk for revision 17986 
 
           | 
        
        
          | 
            File size:
            931 bytes
           | 
        
      
      
| Line |   | 
|---|
| 1 | import numpy as npy
 | 
|---|
| 2 | from TMeltingPoint  import TMeltingPoint
 | 
|---|
| 3 | 
 | 
|---|
| 4 | def DepthAvgTempCond(md):
 | 
|---|
| 5 |    ''' compute conduction dependent temperature profile for an ice sheet. 
 | 
|---|
| 6 |    Usage:
 | 
|---|
| 7 |    Tbar=DepthAvgTempCond(md)
 | 
|---|
| 8 |    '''
 | 
|---|
| 9 | 
 | 
|---|
| 10 |    Tpmp=TMeltingPoint(md.materials.meltingpoint,0) #pressure melting point at 0 pressure.
 | 
|---|
| 11 | 
 | 
|---|
| 12 |    k=md.materials.thermalconductivity
 | 
|---|
| 13 |    G=md.basalforcings.geothermalflux
 | 
|---|
| 14 |    H=md.geometry.thickness[:,0]
 | 
|---|
| 15 |    Ts=md.initialization.temperature
 | 
|---|
| 16 |    alpha=G*H/k
 | 
|---|
| 17 | 
 | 
|---|
| 18 |    Tbar=npy.zeros(md.mesh.numberofvertices,)
 | 
|---|
| 19 | 
 | 
|---|
| 20 |    #find temperature average when we are below melting point: 
 | 
|---|
| 21 |    pos=npy.nonzero( Ts+alpha < Tpmp)
 | 
|---|
| 22 |    if pos:
 | 
|---|
| 23 |            Tbar[pos]=Ts[pos]+alpha[pos]/2 
 | 
|---|
| 24 | 
 | 
|---|
| 25 |    pos=npy.nonzero( Ts+alpha>= Tpmp)
 | 
|---|
| 26 |    if pos:
 | 
|---|
| 27 |            Tbar[pos]=Tpmp+(Tpmp**2-Ts[pos]**2)/2/alpha[pos]+ Tpmp*(Ts[pos]-Tpmp)/alpha[pos]
 | 
|---|
| 28 |    
 | 
|---|
| 29 |    #on ice shelf, easier: 
 | 
|---|
| 30 |    pos=npy.nonzero(md.mask.groundedice_levelset[0]<=0)
 | 
|---|
| 31 |    if pos:
 | 
|---|
| 32 |            Tbar[pos]=(Ts[pos]+Tpmp)/2
 | 
|---|
| 33 | 
 | 
|---|
| 34 |    return Tbar
 | 
|---|
       
      
  Note:
 See   
TracBrowser
 for help on using the repository browser.