Navigation

Operators and Keywords

Function List:

C++ API

Function File: [Ax,Ay] = msh2m_displacement_smoothing(msh,k)

Displace the boundary of a 2D mesh setting a spring with force/length constant k along each edge and enforcing equilibrium.

This function builds matrices containing the resulting (linearized) equation for x and y coordinates of each mesh node. Boundary conditions enforcing the displacement (Dirichlet type problem) or the force (Neumann type) at the boundary must be added to make the system solvable, e.g.: 

           msh = msh2m_structured_mesh(linspace(0,1,10),  linspace(0,1,10),  1,1:4,"left");
          
           dnodes   = msh2m_nodes_on_sides(msh,1:4);
           varnodes = setdiff([1:columns(msh.p)],dnodes);
           xd     = msh.p(1,dnodes)';
           yd     = msh.p(2,dnodes)';
           dx     = dy    = zeros(columns(msh.p),1);
           dxtot  = dytot = -.5*sin(xd.*yd*pi/2);
           Nsteps = 10;
          
           for ii = 1:Nsteps
            dx(dnodes) = dxtot;
            dy(dnodes) = dytot;
            [Ax,Ay] = msh2m_displacement_smoothing(msh,1);
            dx(varnodes) = Ax(varnodes,varnodes) \ ...
                (-Ax(varnodes,dnodes)*dx(dnodes));
            dy(varnodes) = Ay(varnodes,varnodes) \ ...
                (-Ay(varnodes,dnodes)*dy(dnodes));
            msh.p += [ dx'/Nsteps; dy'/Nsteps ] ;
            triplot(msh.t(1:3,:)',msh.p(1,:)',msh.p(2,:)');
            pause(.01)
           endfor

See also: msh2m_jiggle_mesh.

Demonstration 1

The following code

 msh = msh2m_structured_mesh(linspace(0,1,10),
 			    linspace(0,1,10),
 			    1,1:4,"left");
 dnodes   = msh2m_nodes_on_sides(msh,1:4);
 varnodes = setdiff([1:columns(msh.p)],dnodes);
 
 xd = msh.p(1,dnodes)';
 yd = msh.p(2,dnodes)';
 
 dy = zeros(columns(msh.p),1);
 dx = dy;
 
 dxtot = -.5*sin(xd.*yd*pi/2);
 dytot = -.5*sin(xd.*yd*pi/2);
 
 Nsteps = 5;
 for ii=1:Nsteps
   
   dx(dnodes) = dxtot;
   dy(dnodes) = dytot;
   
   [Ax,Ay] = msh2m_displacement_smoothing(msh,1);
   
   dx(varnodes) = Ax(varnodes,varnodes) \ ...
       (-Ax(varnodes,dnodes)*dx(dnodes));
   dy(varnodes) = Ay(varnodes,varnodes) \ ...
       (-Ay(varnodes,dnodes)*dy(dnodes));
 
   msh.p(1,:) += dx'/Nsteps;
   msh.p(2,:) += dy'/Nsteps;
 
   if mod(ii,2)==0
     triplot(msh.t(1:3,:)',msh.p(1,:)',msh.p(2,:)');
     pause(.01)
   endif
 endfor

Produces the following figure

Figure 1

Package: msh