gusto.complex_proxy package
Submodules
gusto.complex_proxy.common module
- gusto.complex_proxy.common.ComplexConstant(z, zi=None)[source]
A namedtuple of 2 firedrake.Constant representing a complex number.
- gusto.complex_proxy.common.ComplexConstantImpl
alias of
ComplexConstant
gusto.complex_proxy.mixed module
gusto.complex_proxy.mixed_impl module
- gusto.complex_proxy.mixed_impl.BilinearForm(W, z, A, return_z=False)[source]
Return a bilinear Form on the complex FunctionSpace W equal to a complex multiple of a bilinear Form on the real FunctionSpace. If z = zr + i*zi is a complex number, u = ur + i*ui is a complex TrialFunction, and b = br + i*bi is a complex linear Form, we want to construct a Form such that (zA)u=b
- (zA)u = (zr*A + i*zi*A)(ur + i*ui)
= (zr*A*ur - zi*A*ui) + i*(zr*A*ui + zi*A*ur)
- = | zr*A -zi*A | | ur | = | br |
- | | | | |zi*A zr*A | | ui | = | bi |
- Parameters:
W – the complex-proxy FunctionSpace
z – a complex number.
A – a generator function for a bilinear Form on the real FunctionSpace, callable as A(*u, *v) where u and v are TrialFunctions and TestFunctions on the real FunctionSpace.
return_z – If true, return Constants for the real/imaginary parts of z used in the BilinearForm.
- gusto.complex_proxy.mixed_impl.ComplexConstant(z, zi=None)[source]
A namedtuple of 2 firedrake.Constant representing a complex number.
- gusto.complex_proxy.mixed_impl.DirichletBC(W, V, bc, function_arg=None)[source]
Return a DirichletBC on the complex FunctionSpace W that is equivalent to the DirichletBC bc on the real FunctionSpace V that W was constructed from.
- Parameters:
W – the complex FunctionSpace.
V – the real FunctionSpace that W was constructed from.
bc – a DirichletBC on the real FunctionSpace that W was constructed from.
- gusto.complex_proxy.mixed_impl.FiniteElement(elem)[source]
Return a UFL FiniteElement which proxies a complex version of the real-valued UFL FiniteElement elem.
The returned complex-valued element has twice as many components as the real-valued element, with each component of the real-valued element having a corresponding ‘real’ and ‘imaginary’ part eg: Non-mixed real elements become 2-component MixedElements. Mixed real elements become MixedElements with 2*len(elem.num_sub_elements()) components. Nested MixedElements are flattened before being proxied.
- Parameters:
elem – the UFL FiniteElement to be proxied
- gusto.complex_proxy.mixed_impl.FunctionSpace(V)[source]
Return a FunctionSpace which proxies a complex version of the real FunctionSpace V.
The returned complex-valued FunctionSpace has twice as many components as the real-valued element, with each component of the real-valued FunctionSpace having a corresponding ‘real’ and ‘imaginary’ part eg: Non-mixed real FunctionSpaces become 2-component MixedFunctionSpaces. Mixed real FunctionSpaces become MixedFunctionSpaces with 2*len(V.ufl_element().num_sub_elements()) components. Function spaces with nested MixedElements are flattened before being proxied.
- Parameters:
V – the real-valued FunctionSpace.
- gusto.complex_proxy.mixed_impl.LinearForm(W, z, f, return_z=False)[source]
Return a Linear Form on the complex FunctionSpace W equal to a complex multiple of a linear Form on the real FunctionSpace. If z = zr + i*zi is a complex number, v = vr + i*vi is a complex TestFunction, we want to construct the Form: <zr*vr,f> + i<zi*vi,f>
- Parameters:
W – the complex-proxy FunctionSpace.
z – a complex number.
f – a generator function for a linear Form on the real FunctionSpace, callable as f(*v) where v are TestFunctions on the real FunctionSpace.
return_z – If true, return Constants for the real/imaginary parts of z used in the LinearForm.
- gusto.complex_proxy.mixed_impl.derivative(z, F, u, return_z=False)[source]
Return a bilinear Form equivalent to z*J where z is a complex number, J = dF/dw, F is a nonlinear Form on the real-valued space, and w is a Function in the real-valued space. The real and imaginary components of the complex Function u most both be equal to w for this operation to be valid.
If z = zr + i*zi is a complex number, x = xr + i*xi is a complex Function, b = br + i*bi is a complex linear Form, J is the bilinear Form dF/dw, we want to construct a Form such that (zJ)x=b
- (zJ)x = (zr*J + i*zi*J)(xr + i*xi)
= (zr*J*xr - zi*J*xi) + i*(zr*A*xi + zi*A*xr)
- = | zr*J -zi*J | | xr | = | br |
- | | | | |zi*J zr*J | | xi | = | bi |
- Parameters:
z – a complex number.
F – a generator function for a nonlinear Form on the real FunctionSpace, callable as F(*u, *v) where u and v are Functions and TestFunctions on the real FunctionSpace.
u – the Function to differentiate F with respect to
return_z – If true, return Constants for the real/imaginary parts of z used in the BilinearForm.
- gusto.complex_proxy.mixed_impl.get_imag(u, vout, name=None)[source]
Copy the imaginary component of the complex Function u into the real-valued Function vout
- Parameters:
u – a complex Function.
vout – A real-valued Function that imaginary component of u is copied into.
- gusto.complex_proxy.mixed_impl.get_real(u, vout)[source]
Copy the real component of the complex Function u into the real-valued Function vout
- Parameters:
u – a complex Function.
vout – A real-valued Function that real component of u is copied into.
- gusto.complex_proxy.mixed_impl.set_imag(u, vnew)[source]
Copy the real-valued Function vnew into the imaginary part of the complex Function u.
- Parameters:
u – a complex Function.
vnew – A real-value Function.
- gusto.complex_proxy.mixed_impl.set_real(u, vnew)[source]
Copy the real-valued Function vnew into the real part of the complex Function u.
- Parameters:
u – a complex Function.
vnew – A real-value Function.
- gusto.complex_proxy.mixed_impl.split(u, i)[source]
- If u is a Coefficient or Argument in the complex FunctionSpace,
returns a tuple with the Function components corresponding to the real or imaginary subelements, indexed appropriately. Analogous to firedrake.split(u)
- Parameters:
u – a Coefficient or Argument in the complex FunctionSpace
i – the index of the components, Part.Real for real or Part.Imag for imaginary.
gusto.complex_proxy.vector module
gusto.complex_proxy.vector_impl module
- gusto.complex_proxy.vector_impl.BilinearForm(W, z, A, return_z=False)[source]
Return a bilinear Form on the complex FunctionSpace W equal to a complex multiple of a bilinear Form on the real FunctionSpace. If z = zr + i*zi is a complex number, u = ur + i*ui is a complex TrialFunction, and b = br + i*bi is a complex linear Form, we want to construct a Form such that (zA)u=b
- (zA)u = (zr*A + i*zi*A)(ur + i*ui)
= (zr*A*ur - zi*A*ui) + i*(zr*A*ui + zi*A*ur)
- = | zr*A -zi*A | | ur | = | br |
- | | | | |zi*A zr*A | | ui | = | bi |
- Parameters:
W – the complex-proxy FunctionSpace
z – a complex number.
A – a generator function for a bilinear Form on the real FunctionSpace, callable as A(*u, *v) where u and v are TrialFunctions and TestFunctions on the real FunctionSpace.
return_z – If true, return Constants for the real/imaginary parts of z used in the BilinearForm.
- gusto.complex_proxy.vector_impl.ComplexConstant(z, zi=None)[source]
A namedtuple of 2 firedrake.Constant representing a complex number.
- gusto.complex_proxy.vector_impl.DirichletBC(W, V, bc, function_arg=None)[source]
Return a DirichletBC on the complex FunctionSpace W that is equivalent to the DirichletBC bc on the real FunctionSpace V that W was constructed from.
- Parameters:
W – the complex FunctionSpace.
V – the real FunctionSpace that W was constructed from.
bc – a DirichletBC on the real FunctionSpace that W was constructed from.
- gusto.complex_proxy.vector_impl.FiniteElement(elem)[source]
Return a UFL FiniteElement which proxies a complex version of the real UFL FiniteElement elem.
The returned complex-valued element has as many components as the real-valued element, but each component has a ‘real’ and ‘imaginary’ part eg: Scalar real elements become 2-vector complex elements. n-vector real elements become 2xn-tensor complex elements (shape)-tensor real elements become (2,shape)-tensor complex elements
- Parameters:
elem – the UFL FiniteElement to be proxied
- gusto.complex_proxy.vector_impl.FunctionSpace(V)[source]
Return a FunctionSpace which proxies a complex version of the real FunctionSpace V.
The returned complex-valued Function space has as many components as the real-valued FunctionSpace, but each component has a ‘real’ and ‘imaginary’ part eg: Scalar components of the real-valued FunctionSpace become 2-vector components of the complex-valued space. n-vector components of the real-valued FunctionSpace become 2xn-tensor components of the complex-valued space. (shape)-tensor components of the real-valued FunctionSpace become (2,shape)-tensor components of the complex-valued FunctionSpace.
- Parameters:
V – the real-valued FunctionSpace.
- gusto.complex_proxy.vector_impl.LinearForm(W, z, f, return_z=False)[source]
Return a Linear Form on the complex FunctionSpace W equal to a complex multiple of a linear Form on the real FunctionSpace. If z = zr + i*zi is a complex number, v = vr + i*vi is a complex TestFunction, we want to construct the Form: <zr*vr,f> + i<zi*vi,f>
- Parameters:
W – the complex-proxy FunctionSpace.
z – a complex number.
f – a generator function for a linear Form on the real FunctionSpace, callable as f(*v) where v are TestFunctions on the real FunctionSpace.
return_z – If true, return Constants for the real/imaginary parts of z used in the LinearForm.
- gusto.complex_proxy.vector_impl.derivative(z, F, u, return_z=False)[source]
Return a bilinear Form equivalent to z*J where z is a complex number, J = dF/dw, F is a nonlinear Form on the real-valued space, and w is a Function in the real-valued space. The real and imaginary components of the complex Function u most both be equal to w for this operation to be valid.
If z = zr + i*zi is a complex number, x = xr + i*xi is a complex Function, b = br + i*bi is a complex linear Form, J is the bilinear Form dF/dw, we want to construct a Form such that (zJ)x=b
- (zJ)x = (zr*J + i*zi*J)(xr + i*xi)
= (zr*J*xr - zi*J*xi) + i*(zr*A*xi + zi*A*xr)
- = | zr*J -zi*J | | xr | = | br |
- | | | | |zi*J zr*J | | xi | = | bi |
- Parameters:
z – a complex number.
F – a generator function for a nonlinear Form on the real FunctionSpace, callable as F(*u, *v) where u and v are Functions and TestFunctions on the real FunctionSpace.
u – the Function to differentiate F with respect to
return_z – If true, return Constants for the real/imaginary parts of z used in the BilinearForm.
- gusto.complex_proxy.vector_impl.get_imag(u, vout, name=None)[source]
Return a real Function equal to the imaginary component of the complex Function u.
- Parameters:
u – a complex Function.
vout – If a real Function then the imaginary component of u is placed here. NotImplementedError(If None then a new Function is returned.)
name – If vout is None, the name of the new Function. Ignored if uout is not none.
- gusto.complex_proxy.vector_impl.get_real(u, vout, name=None)[source]
Return a real Function equal to the real component of the complex Function u.
- Parameters:
u – a complex Function.
vout – If a real Function then real component of u is placed here. NotImplementedError(If None then a new Function is returned.)
name – If vout is None, the name of the new Function. Ignored if vout is not none.
- gusto.complex_proxy.vector_impl.set_imag(u, vnew)[source]
Set the imaginary component of the complex Function u to the value of the real Function vnew.
- Parameters:
u – a complex Function.
vnew – A real Function.
- gusto.complex_proxy.vector_impl.set_real(u, vnew)[source]
Set the real component of the complex Function u to the value of the real Function vnew.
- Parameters:
u – a complex Function.
vnew – A real Function.
- gusto.complex_proxy.vector_impl.split(u, i)[source]
- If u is a Coefficient or Argument in the complex FunctionSpace,
returns a tuple with the Function components corresponding to the real or imaginary subelements, indexed appropriately.
- Parameters:
u – a Coefficient or Argument in the complex FunctionSpace
i – Part.Real for real subelements, Part.Imag for imaginary elements