# [petsc-users] DMPlex with spring elements

Abhyankar, Shrirang G. abhyshr at mcs.anl.gov
Thu Sep 25 12:38:36 CDT 2014

```
From: Miguel Angel Salazar de Troya <salazardetroya at gmail.com<mailto:salazardetroya at gmail.com>>
Date: Thu, 25 Sep 2014 11:43:13 -0500
To: Shri <abhyshr at mcs.anl.gov<mailto:abhyshr at mcs.anl.gov>>
Cc: <petsc-users at mcs.anl.gov<mailto:petsc-users at mcs.anl.gov>>
Subject: Re: [petsc-users] DMPlex with spring elements

I see, and I guess I would have to assign a value of one to the diagonal entry of that degree of freedom in the Jacobian right?

Yes.

Wouldn't this break the symmetry of the Jacobian (in case it were symmetric)?

No.

Shri

Thanks
Miguel

On Sep 25, 2014 11:32 AM, "Abhyankar, Shrirang G." <abhyshr at mcs.anl.gov<mailto:abhyshr at mcs.anl.gov>> wrote:
The solver does not know anything about the boundary conditions. You would have to specify it to the solver by describing the appropriate equations. For e.g. in the power grid example, there is a part in the residual evaluation

if (bus->ide == REF_BUS || bus->ide == ISOLATED_BUS) {
farr[offset] = 0.0;
farr[offset+1] = 0.0;
break;
}

This sets the residual at the nodes marked with REF_BUS or ISOLATED_BUS to 0.0. You can do something similar.

Shri

From: Miguel Angel Salazar de Troya <salazardetroya at gmail.com<mailto:salazardetroya at gmail.com>>
Date: Thu, 25 Sep 2014 10:52:16 -0500
To: Shri <abhyshr at mcs.anl.gov<mailto:abhyshr at mcs.anl.gov>>
Cc: "petsc-users at mcs.anl.gov<mailto:petsc-users at mcs.anl.gov>" <petsc-users at mcs.anl.gov<mailto:petsc-users at mcs.anl.gov>>
Subject: Re: [petsc-users] DMPlex with spring elements

Thanks. Once I have marked the nodes that are fixed nodes using the component data structure, how can I process it later? I mean, at what point does the solver know that those degrees of freedom are actually fixed and how I can tell it that they are fixed?

Miguel

On Thu, Sep 25, 2014 at 10:27 AM, Abhyankar, Shrirang G. <abhyshr at mcs.anl.gov<mailto:abhyshr at mcs.anl.gov>> wrote:

>Thanks. I think the term "Component" was confusing me, I thought it was
>related to the components of a field. I think this would be useful to me
>if I wanted to assign coordinates to the vertices, wouldn't it?

Yes. You can put whatever data you want in the component data structure.

>
>Also, I was wondering how to set up dirichlet boundary conditions,
>basically fixing certain nodes position.
>

>
>
You can add a component at each node with a field marking whether the node
is a boundary node.

>Could I do it as the function SetInitialValues does it in the pflow
>example?
>

No. You need to put in the component data structure before calling

>These values are used to eliminate the zeroth-order energy modes of the
>stiffness matrix?
>

>
>Last question, in my case I have two degrees of freedom per node, when I
>grab the offset with DMNetworkVariableOffset, that's for the first degree
>of freedom in that node and the second degree of freedom would just be
>offset+1?
>

Yes.

Shri

>
>Miguel
>
>
>On Wed, Sep 24, 2014 at 9:52 PM, Abhyankar, Shrirang G.
><abhyshr at mcs.anl.gov<mailto:abhyshr at mcs.anl.gov>> wrote:
>
>If you have equations only at the nodes, with a part of it contributed by
>the edges (springs),  then you can use DMNetwork. If you are planning to
>have equations for the beads in the future, or other higher layers, then
>DMPlex has better functionality
> to manage that.
>
>Shri
>
>
>From: Miguel Angel Salazar de Troya <salazardetroya at gmail.com<mailto:salazardetroya at gmail.com>>
>Date: Wed, 24 Sep 2014 17:38:11 -0500
>To: Shri <abhyshr at mcs.anl.gov<mailto:abhyshr at mcs.anl.gov>>
>Cc: "petsc-users at mcs.anl.gov<mailto:petsc-users at mcs.anl.gov>" <petsc-users at mcs.anl.gov<mailto:petsc-users at mcs.anl.gov>>
>Subject: Re: [petsc-users] DMPlex with spring elements
>
>
>
>
>
>Thanks for your response. I'm attaching a pdf with a description of the
>model. The description of the PetscSection is necessary for the
>DMNetwork? It looks like DMNetwork does not use a PetscSection.
>
>
>
>
>
>
>Miguel
>
>
>On Wed, Sep 24, 2014 at 1:43 PM, Abhyankar, Shrirang G.
><abhyshr at mcs.anl.gov<mailto:abhyshr at mcs.anl.gov>> wrote:
>
>
>>Thanks for your response. My discretization is based on spring elements.
>>For the linear one dimensional case in which each spring has a
>>coefficient k, their jacobian would be this two by two matrix.
>>[  k    -k ]
>>[ -k     k ]
>>
>>and the internal force
>>
>>[ k ( Ui - Uj) ]
>>[ k ( Uj - Ui) ]
>>
>>where Ui and Uj are the node displacements (just one displacement per
>>node because it's one dimensional)
>>
>>For the two dimensional case, assuming small deformations, we have a
>>four-by-four matrix. Each node has two degrees of freedom. We obtain it
>>by performing the outer product of the vector (t , -t) where "t" is the
>>vector that connects both nodes in a spring. This is for the case of
>>small deformations. I would need to assemble each spring contribution to
>>the jacobian and the residual like they were finite elements. The springs
>>share nodes, that's how they are connected. This example is just the
>>linear case, I will have to implement a nonlinear case in a similar
>>fashion.
>>
>>Seeing the DMNetwork example, I think it's what I need, although I don't
>>know much of power electric grids and it's hard for me to understand
>>what's going on. Do you have a good reference to be able to follow the
>>code?
>
>>
>Please see the attached document which has more description of DMNetwork
>and the equations for the power grid example. I don't have anything that
>describes how the power grid example is implemented.
>
>>For example, why are they adding components to the edges?
>>
>><http://www.mcs.anl.gov/petsc/petsc-as/petsc-current/docs/manualpages/DM/
>>D
>>[
>>0],&pfdata.branch[i-eStart]);Miguel
>
>Each edge or node can have several components (limited to 10) attached to
>it. The term components, taken from the circuit terminology, refers to the
>elements of a network. For example, a component could be a resistor,
>inductor, spring, or even edge/vertex weights (for graph problems). For
>code implementation, component is a data structure that holds the data
>needed for the residual, Jacobian, or any other function evaluation. In
>the case of power grid, there are 4 components: branches or transmission
>lines connecting nodes, buses or nodes, generators that are incident at a
>subset of the nodes, and loads that are also incident at a subset of the
>nodes. Each of the these components are defined by their data structures
>given in pf.h.
>
>DMNetwork is a wrapper class of DMPlex specifically for network
>applications that can be solely described using nodes, edges, and their
>associated components. If you have a PDE, or need FEM, or need other
>advanced features then DMPlex would be suitable. Please send us a write-up
>of your equations so that we can assist you better.
>
>Shri
>
>
>>
>>
>>On Tue, Sep 23, 2014 at 11:13 PM, Abhyankar, Shrirang G.
>><abhyshr at mcs.anl.gov<mailto:abhyshr at mcs.anl.gov>> wrote:
>>
>>You may also want to take a look at the DMNetwork framework that can be
>>used for general unstructured networks that don't use PDEs. Its
>>description is given in the manual and an example is in
>>src/snes/examples/tutorials/network/pflow.
>>
>>Shri
>>
>>From:  Matthew Knepley <knepley at gmail.com<mailto:knepley at gmail.com>>
>>Date:  Tue, 23 Sep 2014 22:40:52 -0400
>>To:  Miguel Angel Salazar de Troya <salazardetroya at gmail.com<mailto:salazardetroya at gmail.com>>
>>Cc:  "petsc-users at mcs.anl.gov<mailto:petsc-users at mcs.anl.gov>" <petsc-users at mcs.anl.gov<mailto:petsc-users at mcs.anl.gov>>
>>Subject:  Re: [petsc-users] DMPlex with spring elements
>>
>>
>>>On Tue, Sep 23, 2014 at 4:01 PM, Miguel Angel Salazar de Troya
>>><salazardetroya at gmail.com<mailto:salazardetroya at gmail.com>> wrote:
>>>
>>>Hi all
>>>I was wondering if it could be possible to build a model similar to the
>>>example snes/ex12.c, but with spring elements (for elasticity) instead
>>>of
>>>simplicial elements. Spring elements in a grid, therefore each element
>>>would have two nodes and each node two components. There would be more
>>>differences, because instead of calling the functions f0,f1,g0,g1,g2 and
>>>g3 to build the residual and the jacobian, I would call a routine that
>>>would build the residual vector and the jacobian matrix directly. I
>>>would
>>>not have shape functions whatsoever. My problem is discrete, I don't
>>>have
>>>a PDE and my equations are algebraic. What is the best way in petsc to
>>>solve this problem? Is there any example that I can follow? Thanks in
>>>
>>>
>>>
>>>
>>>Yes, ex12 is fairly specific to FEM. However, I think the right tools
>>>for
>>>what you want are
>>>DMPlex and PetscSection. Here is how I would proceed:
>>>
>>>  1) Make a DMPlex that encodes a simple network that you wish to
>>>simulate
>>>
>>>  2) Make a PetscSection that gets the data layout right. Its hard from
>>>the above
>>>      for me to understand where you degrees of freedom actually are.
>>>This is usually
>>>      the hard part.
>>>
>>>  3) Calculate the residual, so you can check an exact solution. Here
>>>you
>>>use the
>>>      PetscSectionGetDof/Offset() for each mesh piece that you are
>>>interested in. Again,
>>>      its hard to be more specific when I do not understand your
>>>discretization.
>>>
>>>  Thanks,
>>>
>>>    Matt
>>>
>>>
>>>Miguel
>>>
>>>
>>>
>>>--
>>>Miguel Angel Salazar de Troya
>>>Department of Mechanical Science and Engineering
>>>University of Illinois at Urbana-Champaign
>>
>>
>>>(217) 550-2360 <tel:%28217%29%20550-2360>
>>>salaza11 at illinois.edu<mailto:salaza11 at illinois.edu>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>--
>>>What most experimenters take for granted before they begin their
>>>experiments is infinitely more interesting than any results to which
>>>-- Norbert Wiener
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>--
>>Miguel Angel Salazar de Troya
>>Department of Mechanical Science and Engineering
>>University of Illinois at Urbana-Champaign
>>(217) 550-2360
>>salaza11 at illinois.edu<mailto:salaza11 at illinois.edu>
>
>
>
>
>
>
>
>
>
>--
>Miguel Angel Salazar de Troya
>Department of Mechanical Science and Engineering
>University of Illinois at Urbana-Champaign
>(217) 550-2360
>salaza11 at illinois.edu<mailto:salaza11 at illinois.edu>
>
>
>
>
>
>
>
>
>
>
>
>
>
>--
>Miguel Angel Salazar de Troya
>Department of Mechanical Science and Engineering
>University of Illinois at Urbana-Champaign
>(217) 550-2360
>salaza11 at illinois.edu<mailto:salaza11 at illinois.edu>

--
Miguel Angel Salazar de Troya