Re: [firedrake] Varying coefficients in space
Hi Henrik, Let's assume that you have some function mydata(X) which takes as input a npoints * gdim array, where npoints is the number of points at which I need your data values, and gdim is the geometric dimension of the mesh. mydata(X) returns a npoints long vector of the scalar values of mydata evaluated at the points provided. Presumably mydata works by interpolating your external data source, but the precise details are beyond Firedrake's scope. Let's further assume that the Firedrake FunctionSpace into which you wish to interpolate is scalar valued (although this procedure can be extended to vector or tensor valued fields). We'll call that FunctionSpace variable fs. # First, grab the mesh. m = fs.ufl_domain() # Now make the VectorFunctionSpace corresponding to fs. vfs = VectorFunctionSpace(m, fs.ufl_element()) # Next, interpolate the coordinates onto the nodes of vfs. X = interpolate(m.coordinates, vfs) # Make an output function. f = Function(fs) # Use the external data function to interpolate the values of f. f.dat.data[:] = mydata(X.dat.data_ro) That's it. This will also work in parallel, as the interpolation will occur on each process, and Firedrake will take care of the halo updates before the next operation using f. Regards, David On Thu, 6 Apr 2017 at 11:04 Buesing, Henrik <HBuesing@eonerc.rwth-aachen.de> wrote:
Dear all,
I would like to assign a different value for every cell for one of the coefficients of my pde.
I have read all the values for my current discretization with numpy.genfromtxt into a variable. I also have the coefficient Function available.
How can I assign now the values to the Function?
Thank you!
Henrik
--
Dipl.-Math. Henrik Büsing
Institute for Applied Geophysics and Geothermal Energy
E.ON Energy Research Center
RWTH Aachen University
------------------------------------------------------
Mathieustr. 10 | Tel +49 (0)241 80 49907 <+49%20241%208049907>
52074 Aachen, Germany | Fax +49 (0)241 80 49889 <+49%20241%208049889>
------------------------------------------------------
http://www.eonerc.rwth-aachen.de/GGE
hbuesing@eonerc.rwth-aachen.de
------------------------------------------------------
-- Dr David Ham Department of Mathematics Imperial College London
Dear David, Thank you for the explanation! This works! I was missing the dat.data part. Since, I only want to interpolate on DG0, I think it is enough to: # My array which is npoints long perm = numpy.genfromtxt('perm.txt') # My function, constant on each element. K = Function(DG0) # Get the local ranges Kvec = K.vector() lrange = Kvec.local_range() # Assign the values K.dat.data[:] = perm[lrange[0]:lrange[1]] Thank you! Henrik -- Dipl.-Math. Henrik Büsing Institute for Applied Geophysics and Geothermal Energy E.ON Energy Research Center RWTH Aachen University ------------------------------------------------------ Mathieustr. 10 | Tel +49 (0)241 80 49907 52074 Aachen, Germany | Fax +49 (0)241 80 49889 ------------------------------------------------------ http://www.eonerc.rwth-aachen.de/GGE hbuesing@eonerc.rwth-aachen.de<mailto:hbuesing@eonerc.rwth-aachen.de> ------------------------------------------------------ Von: firedrake-bounces@imperial.ac.uk<mailto:firedrake-bounces@imperial.ac.uk> [mailto:firedrake-bounces@imperial.ac.uk] Im Auftrag von David Ham Gesendet: 06 April 2017 12:38 An: firedrake <firedrake@imperial.ac.uk<mailto:firedrake@imperial.ac.uk>> Betreff: Re: [firedrake] Varying coefficients in space Hi Henrik, Let's assume that you have some function mydata(X) which takes as input a npoints * gdim array, where npoints is the number of points at which I need your data values, and gdim is the geometric dimension of the mesh. mydata(X) returns a npoints long vector of the scalar values of mydata evaluated at the points provided. Presumably mydata works by interpolating your external data source, but the precise details are beyond Firedrake's scope. Let's further assume that the Firedrake FunctionSpace into which you wish to interpolate is scalar valued (although this procedure can be extended to vector or tensor valued fields). We'll call that FunctionSpace variable fs. # First, grab the mesh. m = fs.ufl_domain() # Now make the VectorFunctionSpace corresponding to fs. vfs = VectorFunctionSpace(m, fs.ufl_element()) # Next, interpolate the coordinates onto the nodes of vfs. X = interpolate(m.coordinates, vfs) # Make an output function. f = Function(fs) # Use the external data function to interpolate the values of f. f.dat.data[:] = mydata(X.dat.data_ro) That's it. This will also work in parallel, as the interpolation will occur on each process, and Firedrake will take care of the halo updates before the next operation using f. Regards, David On Thu, 6 Apr 2017 at 11:04 Buesing, Henrik <HBuesing@eonerc.rwth-aachen.de<mailto:HBuesing@eonerc.rwth-aachen.de>> wrote: Dear all, I would like to assign a different value for every cell for one of the coefficients of my pde. I have read all the values for my current discretization with numpy.genfromtxt into a variable. I also have the coefficient Function available. How can I assign now the values to the Function? Thank you! Henrik -- Dipl.-Math. Henrik Büsing Institute for Applied Geophysics and Geothermal Energy E.ON Energy Research Center RWTH Aachen University ------------------------------------------------------ Mathieustr. 10 | Tel +49 (0)241 80 49907<tel:+49%20241%208049907> 52074 Aachen, Germany | Fax +49 (0)241 80 49889<tel:+49%20241%208049889> ------------------------------------------------------ http://www.eonerc.rwth-aachen.de/GGE hbuesing@eonerc.rwth-aachen.de<mailto:hbuesing@eonerc.rwth-aachen.de> ------------------------------------------------------ -- Dr David Ham Department of Mathematics Imperial College London
On 06/04/17 15:49, Buesing, Henrik wrote:
Dear David,
Thank you for the explanation! This works! I was missing the dat.data part. Since, I only want to interpolate on DG0, I think it is enough to:
# My array which is npoints long
perm = numpy.genfromtxt('perm.txt')
# My function, constant on each element. K = Function(DG0)
# Get the local ranges Kvec = K.vector()
lrange = Kvec.local_range()
# Assign the values
K.dat.data[:] = perm[lrange[0]:lrange[1]]
This assumes that the data you've made are in the same ordering as the function space degree of freedom numbering. We make no guarantees not to change this numbering without warning behind your back, so I recommend not using this approach. Lawrence
participants (3)
-
Buesing, Henrik
-
David Ham
-
Lawrence Mitchell