Where can i find the dolfin module files?

General
1

Hello everyone,

I installed fenics Legacy with Docker (windows 11) following this tutorial :
https://fenics.readthedocs.io/projects/containers/en/latest/jupyter.html

i runned this command so i can use fenics with jupyter:

docker run --name notebook -w /home/fenics -v ${pwd}:/home/fenics/shared -d -p 127.0.0.1:8888:8888 quay.io/fenicsproject/stable 'jupyter-notebook --ip=0.0.0.0'

I would like to open the files from the dolfin module, especially the dolfin/io/XDMFFile.cpp file… I looked pretty much everywhere in the Docker container’s files but can’t find the directory containing these files.
I am an absolute begginer in using Docker and fenics, any advice on where can i find these files would be truly appreciated !

2

The legacy images remove the source code after installation, Ref: https://bitbucket.org/fenics-project/docker/raw/b1d5e0963afc9a23d11a18b1e82347cee21ae520/dockerfiles/stable/Dockerfile
This is Also done in the more modern images:
Package fenics-gmsh · GitHub
However, it is quite easy to take one of the modern docker files, and remove a single line at the end of the file: https://github.com/scientificcomputing/packages/blob/main/fenics/Dockerfile#L202 and build the image yourself.

You could simply clone dolfin from
Bitbucket and inspect the files locally. What exactly are you looking for?

3

First of all, thanks you for your reply.

I will precise what i am trying to achieve.
I have the exact same issue covered in this topic to which you responded :
https://fenicsproject.discourse.group/t/cannot-open-xdmf-files-in-legacy-version-of-fenics-latest-paraview/10344

My understanding of your last answer is that i have to modify the 5 files as it is presented in this link https://bitbucket.org/fenics-project/dolfin/commits/6fbc9fb8e1272c4cdd73d7e2e432dc29ebfe03ac
That’s why i’m searching those files in order to edit them.
I’m totally new on these topics so please correct me if my understanding is wrong.

4

This has been fixed in fenics versions · scientificcomputing · GitHub
or if you want jupyter-lab support:
fenics-gmsh versions · scientificcomputing · GitHub
were you run

docker run -ti -v $(pwd):/root/shared -w /root/shared --name=notebook -p 8888:8888 ghcr.io/scientificcomputing/fenics-gmsh:jupterlab-2023-04-21 

jupyter lab --ip 0.0.0.0 --no-browser --allow-root
5

Thanks again for your reply. I tried but i still have the same issue, here is what i did:
I runned the windows terminal as admin and runned the two commands :

docker run -ti -v $(pwd):/root/shared -w /root/shared --name=notebook -p 8888:8888 ghcr.io/scientificcomputing/fenics-gmsh:jupterlab-2023-04-21 

jupyter lab --ip 0.0.0.0 --no-browser --allow-root

Everything worked well. Then i opened a new jupyter file and pasted this code :

import fenics as fe

CANTILEVER_LENGTH = 1.0
CANTILEVER_WIDTH = 0.2

N_POINTS_LENGTH = 10
N_POINTS_WIDTH = 3

LAME_MU = 1.0
LAME_LAMBDA = 1.25
DENSITY = 1.0
ACCELERATION_DUE_TO_GRAVITY = 0.016

def main():
    # Mesh and Vector Function Space
    mesh = fe.BoxMesh(
        fe.Point(0.0, 0.0, 0.0),
        fe.Point(CANTILEVER_LENGTH, CANTILEVER_WIDTH, CANTILEVER_WIDTH),
        N_POINTS_LENGTH,
        N_POINTS_WIDTH,
        N_POINTS_WIDTH,
    )
    lagrange_vector_space_first_order = fe.VectorFunctionSpace(
        mesh,
        "Lagrange",
        1,
    )
    
    # Boundary Conditions
    def clamped_boundary(x, on_boundary):
        return on_boundary and x[0] < fe.DOLFIN_EPS

    dirichlet_clamped_boundary = fe.DirichletBC(
        lagrange_vector_space_first_order,
        fe.Constant((0.0, 0.0, 0.0)),
        clamped_boundary,
    )

    # Define strain and stress
    def epsilon(u):
        engineering_strain = 0.5 * (fe.nabla_grad(u) + fe.nabla_grad(u).T)
        return engineering_strain
    
    def sigma(u):
        cauchy_stress = (
            LAME_LAMBDA * fe.tr(epsilon(u)) * fe.Identity(3)
            +
            2 * LAME_MU * epsilon(u)
        )
        return cauchy_stress
    
    # Define weak form
    u_trial = fe.TrialFunction(lagrange_vector_space_first_order)
    v_test = fe.TestFunction(lagrange_vector_space_first_order)
    forcing = fe.Constant((0.0, 0.0, - DENSITY * ACCELERATION_DUE_TO_GRAVITY))
    traction = fe.Constant((0.0, 0.0, 0.0))

    weak_form_lhs = fe.inner(sigma(u_trial), epsilon(v_test)) * fe.dx  # Crucial to use inner and not dot
    weak_form_rhs = (
        fe.dot(forcing, v_test) * fe.dx
        +
        fe.dot(traction, v_test) * fe.ds
    )

    # Compute solution
    u_solution = fe.Function(lagrange_vector_space_first_order)
    fe.solve(
        weak_form_lhs == weak_form_rhs,
        u_solution,
        dirichlet_clamped_boundary,
    )

    # Compute the von Mises stress
    deviatoric_stress_tensor = (
        sigma(u_solution)
        -
        1/3 * fe.tr(sigma(u_solution)) * fe.Identity(3)
    )
    von_Mises_stress = fe.sqrt(3/2 * fe.inner(deviatoric_stress_tensor, deviatoric_stress_tensor))

    lagrange_scalar_space_first_order = fe.FunctionSpace(
        mesh,
        "Lagrange",
        1,
    )
    von_Mises_stress = fe.project(von_Mises_stress, lagrange_scalar_space_first_order)

    # Write out fields for visualization with Paraview
    u_solution.rename("Displacement Vector", "")
    von_Mises_stress.rename("von Mises stress", "")

    beam_deflection_file = fe.XDMFFile("beam_deflection.xdmf")
    beam_deflection_file.parameters["flush_output"] = True
    beam_deflection_file.parameters["functions_share_mesh"] = True
    beam_deflection_file.write(u_solution, 0.0)
    beam_deflection_file.write(von_Mises_stress, 0.0)


if __name__ == "__main__":
    main()

The code compiled without errors.
I tried to opened the .xdmf file that the code created in paraview :

<?xml version="1.0"?>
<!DOCTYPE Xdmf SYSTEM "Xdmf.dtd" []>
<Xdmf Version="3.0" xmlns:xi="http://www.w3.org/2001/XInclude">
  <Domain>
    <Grid Name="TimeSeries" GridType="Collection" CollectionType="Temporal">
      <Grid Name="mesh" GridType="Uniform">
        <Topology NumberOfElements="540" TopologyType="Tetrahedron" NodesPerElement="4">
          <DataItem Dimensions="540 4" NumberType="UInt" Format="HDF">beam_deflection.h5:/Mesh/0/mesh/topology</DataItem>
        </Topology>
        <Geometry GeometryType="XYZ">
          <DataItem Dimensions="176 3" Format="HDF">beam_deflection.h5:/Mesh/0/mesh/geometry</DataItem>
        </Geometry>
        <Time Value="0.000000000000000e+00" />
        <Attribute Name="Displacement Vector" AttributeType="Vector" Center="Node">
          <DataItem Dimensions="176 3" Format="HDF">beam_deflection.h5:/VisualisationVector/0</DataItem>
        </Attribute>
        <Attribute Name="von Mises stress" AttributeType="Scalar" Center="Node">
          <DataItem Dimensions="176 1" Format="HDF">beam_deflection.h5:/VisualisationVector/1</DataItem>
        </Attribute>
      </Grid>
    </Grid>
  </Domain>
</Xdmf>

And paraview crashed when i clicked on apply with both xdmf3ReaderS xdmf3ReaderT. With XDMF Reader, it didn’t crash but i had this error :

ERROR: In vtkXdmfReader.cxx, line 487
vtkXdmfReader (000001E8DABCCA00): Failed to read data.

ERROR: In vtkExecutive.cxx, line 741
vtkPVCompositeDataPipeline (000001E8DAB4C770): Algorithm vtkFileSeriesReader (000001E8DAA2C2B0) returned failure for request: vtkInformation (000001E8E9E9D9D0)
Debug: Off
Modified Time: 411165
Reference Count: 1
Registered Events: (none)
Request: REQUEST_DATA
FORWARD_DIRECTION: 0
ALGORITHM_AFTER_FORWARD: 1
FROM_OUTPUT_PORT: 0

6

The error is at least not with the xdmf file, as the paths are complete

Im not at a laptop, so i cant test the code right now

7

Something weird is that if i add these lines to the code as you suggested in the other topic, the issue disepears :

    encoding = fe.XDMFFile.Encoding.ASCII
    file = fe.XDMFFile("beam_deflection.xdmf")
    file.write(u_solution,encoding=encoding)