In [ ]:

    

# Rayleigh formula for free thermal convection, formulae taken from Kuhn and Gessner 2009, Surveys in Geophysics

def Rayleigh():
    """ Properties of pure water at 20 degrees K """
    # values of pure water taken from http://people.ucsc.edu/~bkdaniel/WaterProperties.html
    alpha =((0.207 + 0.385 + 0.523) / 3) * 1*10**-3 # thermal expansivity of water
    print alpha
    rho = 998.0 # Density of pure water
    c = 4179.0 # Fluid heat capacity
    g = 9.812
    mu = ((1.004 + 0.658 + 0.475) / 3) * 1*10**-3 #dynamic viscosity of water
    """ Mean Values """    
    WLFM0_xyz = 2.7, #np.array(S1_out.get_array_as_xyz_structure("WLFM0"))
    WLFM0_mean = 2.7 #np.mean(WLFM0_xyz[:, :, :]) # Mean of medium's thermal expansivity
    #print WLFM0_mean
    
    perm_xyz = 1E-12 # np.array(S1_out.get_array_as_xyz_structure("PERM"))
    mean_perm = 1E-12 # np.mean(perm_xyz[:,:,:]) # Mean of permeability field
    #mean_perm = 5.835e-12
    #print "%0.4g" %mean_perm
    
   #temp_xyz = np.array(S1_out.get_array_as_xyz_structure("TEMP")) # determine temperature gradient at the subplot
    T2 = 170. #temp_xyz[0:, :, 0]
    T1 = 20. # temp_xyz[0:, :, -1]
    # z_axis = temp_xyz[0, :, 0:]
    # z_axis_list = list(z_axis[:,:].ravel(order = "F"))
    length = 1000. # len(z_axis_list) * 100.0
    #print length
    Value = (T2 - T1) / length
    dtemp_value = Value # np.mean(Value[:, :])
    
    print dtemp_value
    
    Ra = ((T2 - T1) * mean_perm * length * alpha * (rho**2) * g * c) / (WLFM0_mean * mu)
    return Ra
print "%0.4g" %Rayleigh()


    

In [7]:

    

from vtk import *

import matplotlib
matplotlib.use('Agg')
from matplotlib.figure import Figure
from matplotlib.backends.backend_agg import FigureCanvasAgg
import pylab as p

# The vtkImageImporter will treat a python string as a void pointer
importer = vtkImageImport()
importer.SetDataScalarTypeToUnsignedChar()
importer.SetNumberOfScalarComponents(4)

# It's upside-down when loaded, so add a flip filter
imflip = vtkImageFlip()
imflip.SetInputData(importer.GetOutput())
imflip.SetFilteredAxis(1)

# Map the plot as a texture on a cube
cube = vtkCubeSource()

cubeMapper = vtkPolyDataMapper()
cubeMapper.SetInputData(cube.GetOutput())

cubeActor = vtkActor()
cubeActor.SetMapper(cubeMapper)

# Create a texture based off of the image
cubeTexture = vtkTexture()
cubeTexture.InterpolateOn()
cubeTexture.SetInputData(imflip.GetOutput())
cubeActor.SetTexture(cubeTexture)

ren = vtkRenderer()
ren.AddActor(cubeActor)

renWin = vtkRenderWindow()
renWin.AddRenderer(ren)

iren = vtkRenderWindowInteractor()
iren.SetRenderWindow(renWin)

# Now create our plot
fig = Figure()
canvas = FigureCanvasAgg(fig)
ax = fig.add_subplot(111)
ax.grid(True)
ax.set_xlabel('Hello from VTK!', size=16)
ax.bar(xrange(10), p.rand(10))

# Powers of 2 image to be clean
w,h = 1024, 1024
dpi = canvas.figure.get_dpi()
#fig.set_figsize_inches(w / dpi, h / dpi)
canvas.draw() # force a draw

# This is where we tell the image importer about the mpl image
extent = (0, w - 1, 0, h - 1, 0, 0)
importer.SetWholeExtent(extent)
importer.SetDataExtent(extent)
# importer.SetImportVoidPointer(canvas.buffer_rgba(0,0), 1)
importer.Update()

iren.Initialize()
iren.Start()


    

In [ ]:

    

imflip.SetInputData(importer.GetOutput())


    

In [ ]: