In [1]:

    

#Henry Steven Rueda Corredor Tarea 3 Geodinámica
%pylab inline
import math
import numpy as np
import scipy.special as sp
import matplotlib.pyplot as plt


    

    




Populating the interactive namespace from numpy and matplotlib

In [2]:

    

#const   # T1-T0=Ts=1300 Kelvin  
k=10**-6
Ts=-1300
vc=1/(100.0*365.0*24*3600)
densmantle=3300
denswater=1000
alpha=1.6*10**-5
T2=1200
densidad=(densmantle/(densmantle-denswater))
R=8.314
ep=523*10**3
H=4.2*10**5
E=7*10**4


    

In [3]:

    

def Equation(x,y,v):
    T= sp.erfc((abs(y*1000)*(sqrt(v*vc)))/(2*(sqrt(k*abs(x*1000)))))*-Ts+1573
    return T


    

In [4]:

    

x=linspace(-100,100,1000)
y=linspace(0,-60,1000)
X,Y=np.meshgrid(x,y)


    

In [5]:

    

a=Equation (X,Y,0.5)
b=Equation (X,Y,4.0)
c=Equation (X,Y,10.0)


    

In [45]:

    

fig,axes=plt.subplots(1,1,figsize=(20,5))
plt.imshow(a,extent=[x.min(), x.max(), y.min(), y.max()])
plt.title('Temperature for 0.5 cm/year',fontsize=20)
plt.xlabel('x (m)', fontsize=20)
plt.ylabel('y (m)', fontsize=15)
cb=plt.colorbar()


    

    




C:\Users\henry\Anaconda3\lib\site-packages\matplotlib\collections.py:590: FutureWarning: elementwise comparison failed; returning scalar instead, but in the future will perform elementwise comparison
  if self._edgecolors == str('face'):






    

In [44]:

    

fig,axes=plt.subplots(1,1,figsize=(20,5))
plt.imshow(b,extent=[x.min(),x.max(),y.min(),y.max()])
plt.title('Temperature for 4 cm/year',fontsize=20)
plt.xlabel('x (m)', fontsize=20)
plt.ylabel('y (m)', fontsize=15)
cb=plt.colorbar()


    

    




C:\Users\henry\Anaconda3\lib\site-packages\matplotlib\collections.py:590: FutureWarning: elementwise comparison failed; returning scalar instead, but in the future will perform elementwise comparison
  if self._edgecolors == str('face'):






    

In [43]:

    

fig,axes=plt.subplots(1,1,figsize=(20,5))
plt.imshow(c,extent=[x.min(),x.max(),y.min(),y.max()])
plt.title('Temperature for 10 cm/year',fontsize=20)
plt.xlabel('x (m)', fontsize=20)
plt.ylabel('y (m)', fontsize=15)
cb=plt.colorbar()


    

    




C:\Users\henry\Anaconda3\lib\site-packages\matplotlib\collections.py:590: FutureWarning: elementwise comparison failed; returning scalar instead, but in the future will perform elementwise comparison
  if self._edgecolors == str('face'):






    

In [46]:

    

def Thickness(v,x):
    Th= -2.32*(sqrt(k*abs(x*1000)))/(sqrt(v*vc))
    return Th


    

In [47]:

    

d=Thickness(0.5,x)
e=Thickness(4,x)
f=Thickness(10,x)


    

In [48]:

    

plot(x,d)
plot(x,e)
plot(x,f)
title("Thickness",fontsize=15)
xlabel('x (m)', fontsize=20)
ylabel('y (m)', fontsize=15)


    

    
Out[48]:





<matplotlib.text.Text at 0x2c7bd68>






    

In [49]:

    

def Bathymetry(v,x):
    bath= -densidad*alpha*T2*2*(sqrt(k*abs(x*1000)))/(sqrt(v*vc*pi))
    return bath


    

In [50]:

    

g=Bathymetry(0.5,x)
h=Bathymetry(4,x)
i=Bathymetry(10,x)


    

In [51]:

    

plot(x,g)
plot(x,h)
plot(x,i)
title("Bathymetry",fontsize= 10)
xlabel('x(m)', fontsize=15)
ylabel('y(m)',fontsize=15)


    

    
Out[51]:





<matplotlib.text.Text at 0x2b3ff98>






    

In [52]:

    

def Trelax(t,sig):
    Phi=((3*exp(ep/(R*t)))/(2*E*H*sig**2))*(100.0*vc)
    return Phi


    

In [53]:

    

t=linspace(600,2000,500)
na=Trelax(t,10)
ne=Trelax(t,100)
ni=Trelax(t,1000)


    

In [54]:

    

plot(t,na)
plot(t,ne)
plot(t,ni)
plt.yscale('log')
title("Relaxation time", fontsize=15)
xlabel('T(K)')
ylabel('trelax (years)')
xlim(600,2000)
ylim(0,10**10)


    

    
Out[54]:





(7.3924606212830455e-11, 10000000000)






    

In [ ]:

    




    

In [ ]: