### Using "project" to compute the magnetic field from the potential field

341

views

0

Hello everybody,

I have a problem with the differentiation of the potential vector Ato compute the magnetic field B.

I have made the tutorial "magnetostatics" to compute the potential vector generated by straight wires, in this tutorial the instruction:

B = project(as_vector((A_z.dx(1), -A_z.dx(0))), W)

is used to compute B components as the partial derivatives of A.

I have then tried to solve my own simple magnetostatics problem, namely compute the distribution of the magnetic field generated by a single wire like the (Biot Savart law).

THe potential vector A is computed but when I try to get the field B, it seems to be computed but none plot appears and even in ParaView nothing appears.

Which error have I made with "dx", "project" ?

HEre below is my own script

from __future__ import division

from fenics import *

import numpy as np

import matplotlib.pyplot as plt

from mshr import *

# Problem data

I = 100 # [A] current intensity

r = 0.01 # [m] wire radius

A = np.pi*pow(r,2) # [m2] wire cross sectional area

mo = 4*np.pi*1e-7 # [N/A2] air magnetic permeability

mu_r = 1 # [-] relative magnetic permeability

mu = mo*mu_r

# FE solver data

R = 1 # [m] mesh radius - mesh boundary

N_elem = 60 # [-] number of elements created from radius mesh

domain = Circle(Point(0,0),R)

mesh = generate_mesh(domain,100)

J = Constant(I/A) # [A/m2] current density in the wire

V = FunctionSpace(mesh,'P',1)

Az = TrialFunction(V)

v = TestFunction(V)

def boundary(x,on_boundary):

return on_boundary

Az_D = Constant(0)

bc = DirichletBC(V,Az_D,boundary)

dx = Measure('dx',domain=mesh)

a = (1/mu)*dot(grad(Az),grad(v))*dx

L = J*v*dx

Az = Function(V)

solve(a==L,Az,bc)

W = VectorFunctionSpace(mesh,'P',1)

B = project(as_vector((Az.dx(1),-Az.dx(0))),W)

plot(mesh)

interactive()

plot(Az)

interactive()

plot(B)

interactive()

vtkfile = File('Spira_B_solution.pvd')

vtkfile << B

vtkfile_Az = File('Spira_Az_solution.pvd')

vtkfile_Az << Az

Thank you in advance for your help

Diego

I have a problem with the differentiation of the potential vector Ato compute the magnetic field B.

I have made the tutorial "magnetostatics" to compute the potential vector generated by straight wires, in this tutorial the instruction:

B = project(as_vector((A_z.dx(1), -A_z.dx(0))), W)

is used to compute B components as the partial derivatives of A.

I have then tried to solve my own simple magnetostatics problem, namely compute the distribution of the magnetic field generated by a single wire like the (Biot Savart law).

THe potential vector A is computed but when I try to get the field B, it seems to be computed but none plot appears and even in ParaView nothing appears.

Which error have I made with "dx", "project" ?

HEre below is my own script

from __future__ import division

from fenics import *

import numpy as np

import matplotlib.pyplot as plt

from mshr import *

# Problem data

I = 100 # [A] current intensity

r = 0.01 # [m] wire radius

A = np.pi*pow(r,2) # [m2] wire cross sectional area

mo = 4*np.pi*1e-7 # [N/A2] air magnetic permeability

mu_r = 1 # [-] relative magnetic permeability

mu = mo*mu_r

# FE solver data

R = 1 # [m] mesh radius - mesh boundary

N_elem = 60 # [-] number of elements created from radius mesh

domain = Circle(Point(0,0),R)

mesh = generate_mesh(domain,100)

J = Constant(I/A) # [A/m2] current density in the wire

V = FunctionSpace(mesh,'P',1)

Az = TrialFunction(V)

v = TestFunction(V)

def boundary(x,on_boundary):

return on_boundary

Az_D = Constant(0)

bc = DirichletBC(V,Az_D,boundary)

dx = Measure('dx',domain=mesh)

a = (1/mu)*dot(grad(Az),grad(v))*dx

L = J*v*dx

Az = Function(V)

solve(a==L,Az,bc)

W = VectorFunctionSpace(mesh,'P',1)

B = project(as_vector((Az.dx(1),-Az.dx(0))),W)

plot(mesh)

interactive()

plot(Az)

interactive()

plot(B)

interactive()

vtkfile = File('Spira_B_solution.pvd')

vtkfile << B

vtkfile_Az = File('Spira_Az_solution.pvd')

vtkfile_Az << Az

Thank you in advance for your help

Diego

Community: FEniCS Project

### 1 Answer

0

I am sorry I was wrong.

The problem has been solved.

Diego
Hi Diego, would you like to point out how you solved this? Might be helpful to me right now. Thanks in advance.

The problem has been solved.

Diego

written
5 weeks ago by
Rebecca Pflanze

Please login to add an answer/comment or follow this question.