Demo entry 2810420

saasas

   

Submitted by anonymous on Oct 03, 2015 at 15:29
Language: Python. Code size: 3.4 kB.

from dolfin import*
import scipy.special

p = 3

class TFactor01(Expression):
    def eval(self, value, x):
	r = (x[0] ** 2 + x[1] ** 2 + x[2] ** 2) ** 0.5
	rho = (x[0] ** 2 + x[1] ** 2) ** 0.5
	if r < 0.01:
            sinTheta = 0
            cosTheta = 0
        else:
            sinTheta = rho / r
            cosTheta = x[2] / r
	if rho < 0.01:
            sinPhi = 0
            cosPhi = 0
        else:
            sinPhi = x[1] / rho
            cosPhi = x[0] / rho
	eR = [sinTheta * cosPhi, sinTheta * sinPhi, cosTheta]
	t01 = r ** 2 * (1.0 - r ** 2)
	angleFactorT01 = scipy.special.lpmv(0,1,cosTheta)*1.0
	totalFactorT01 = t01 * angleFactorT01
	value[0] = totalFactorT01 * eR[0]
	value[1] = totalFactorT01 * eR[1] 
	value[2] = totalFactorT01 * eR[2]
    def value_shape(self):
	return (3,)
        
class SFactor02(Expression):
    def eval(self, value, x):
	r = (x[0] ** 2 + x[1] ** 2 + x[2] ** 2) ** 0.5
	rho = (x[0] ** 2 + x[1] ** 2) ** 0.5
	if r < 0.01:
            sinTheta = 0
            cosTheta = 0
        else:
            sinTheta = rho / r
            cosTheta = x[2] / r
	if rho < 0.01:
            sinPhi = 0
            cosPhi = 0
        else:
            sinPhi = x[1] / rho
            cosPhi = x[0] / rho
	eR = [sinTheta * cosPhi, sinTheta * sinPhi, cosTheta]
        s02 = r ** 6 * (1.0 - r ** 2) ** 3
        angleFactorS02 = scipy.special.lpmv(0,2,cosTheta)*1.0
        totalFactorS02 = s02 * angleFactorS02
	value[0] = totalFactorS02 * eR[0]
	value[1] = totalFactorS02 * eR[1] 
	value[2] = totalFactorS02 * eR[2]
    def value_shape(self):
	return (3,)
      
class SFactor2c2(Expression):
    def eval(self, value, x):
	r = (x[0] ** 2 + x[1] ** 2 + x[2] ** 2) ** 0.5
	rho = (x[0] ** 2 + x[1] ** 2) ** 0.5
	if r < 0.01:
            sinTheta = 0
            cosTheta = 0
        else:
            sinTheta = rho / r
            cosTheta = x[2] / r
	if rho < 0.01:
            sinPhi = 0
            cosPhi = 0
        else:
            sinPhi = x[1] / rho
            cosPhi = x[0] / rho
	eR = [sinTheta * cosPhi, sinTheta * sinPhi, cosTheta]
        s2c2 = r ** 4 * (1.0 - r ** 2) ** 2 * cos(p*pi*r)
        angleFactorS2c2 = scipy.special.lpmv(2,2,cosTheta)*(cosPhi*cosPhi - sinPhi*sinPhi)
        totalFactorS2c2 = s2c2 * angleFactorS2c2
	value[0] = totalFactorS2c2 * eR[0]
	value[1] = totalFactorS2c2 * eR[1] 
	value[2] = totalFactorS2c2 * eR[2]
    def value_shape(self):
	return (3,)
        
class SFactor2s2(Expression):
    def eval(self, value, x):
	r = (x[0] ** 2 + x[1] ** 2 + x[2] ** 2) ** 0.5
	rho = (x[0] ** 2 + x[1] ** 2) ** 0.5
	if r < 0.01:
            sinTheta = 0
            cosTheta = 0
        else:
            sinTheta = rho / r
            cosTheta = x[2] / r
	if rho < 0.01:
            sinPhi = 0
            cosPhi = 0
        else:
            sinPhi = x[1] / rho
            cosPhi = x[0] / rho
	eR = [sinTheta * cosPhi, sinTheta * sinPhi, cosTheta]
        s2s2 = r ** 4 * (1.0 - r ** 2) ** 2 * sin(p*pi*r)
        angleFactorS2s2 = scipy.special.lpmv(2,2,cosTheta)*(2*sinPhi*cosPhi)
        totalFactorS2s2 = s2s2 * angleFactorS2s2
	value[0] = totalFactorS2s2 * eR[0]
	value[1] = totalFactorS2s2 * eR[1] 
	value[2] = totalFactorS2s2 * eR[2]
    def value_shape(self):
	return (3,)

This snippet took 0.01 seconds to highlight.

Back to the Entry List or Home.

Delete this entry (admin only).