Download :download:`this testcase <2d_arc_boutement.py>`.

Arc-boutement
=============
This example demonstrates the arc-boutement phenomenon in a simple
two-dimensional granular configuration. A small chain of circular particles
is compressed horizontally against a fixed particle, forming a self-locking
structure due to frictional contacts. This effect illustrates how granular
materials can support compressive loads through contact force chain.

Keywords
--------
DEM, Rigid Body, Friction

.. code-block:: python
  
  import sys, os, shutil
  import numpy as np
  from migflow import scontact, time_integration
  

Output Directory
----------------
The output directory is prepared for storing the simulation results.

.. code-block:: python
  
  outputdir = "output" if len(sys.argv) < 2 else sys.argv[1]
  shutil.rmtree(outputdir, ignore_errors=True)
  os.makedirs(outputdir)
  
  

Particle Problem
----------------
We begin by defining the 2D particle problem and configuring the contact
solver. The contact geometry is updated dynamically, and the prediction
basis is activated to improve convergence of the contact iterations.

.. code-block:: python
  
  p = scontact.ParticleProblem(2)
  p.set_fixed_contact_geometry(0)
  p.set_predict_basis(1)
  

Material Properties and Parameters
----------------------------------
The particle density, radius, gravity vector and the friction coefficient are defined.
A large friction coefficient is assigned between particles to ensure interlocking and prevent sliding.

.. code-block:: python
  
  
  rho = 1000
  r = 0.01
  g = np.array([0.0, -9.81])
  friction_coefficient = 5
  p.set_friction_coefficient(friction_coefficient, "particle", "particle")
  p.set_friction_coefficient(0.0, "particle", "wall")

Body Definition
---------------
A small chain of five particles is created and initialized as a rigid body.
The body is placed slightly above the ground, with each particle in contact
with its neighbors. The high inverse inertia value ensures the body behaves
as a nearly rigid object.

.. code-block:: python
  
  n = 5  # number of particles in the chain
  xp = np.array([0, 3 * r])  # body center position
  x_rel = np.arange(n) * 2 * r  # relative particle spacing
  x_rel = np.vstack((x_rel, np.zeros_like(x_rel))).T  # particle coordinates
  radii = np.full(x_rel.shape[0], r)  # uniform radii
  masses = np.full(x_rel.shape[0], np.pi * rho * r**2)  # uniform masses
  body = p.init_body(xp, x_rel, masses, radii, material="particle")
  

Fixed Particle
--------------
A fixed particle is added to act as a support against which the body
will push. This setup allows the arc-boutement effect to develop between
the free body and the fixed particle through frictional contact.

.. code-block:: python
  
  fixed_body = p.add_body(np.array([0.0, r]), 0.0, 0.0)
  p.add_particle_to_body(np.array([0.0, 0.0]), r, fixed_body, material="particle")
  

Boundary Conditions
-------------------
A horizontal wall is added below the particles to prevent vertical motion.

.. code-block:: python
  
  p.add_boundary_segment(
      np.array([-10 * r, 0.0]), np.array([50 * r, 0.0]), tag="wall", material="wall"
  )
  

Time Integration
----------------
The simulation runs for a short time interval with a fixed time step.
At each iteration, the gravitational forces are applied to the particles,
and the system evolves according to the contact dynamics solver.

.. code-block:: python
  
  t = 0.0
  i = 0
  dt = 1e-3
  tEnd = 250 * dt
  outf = 10
  
  while t < tEnd:
      print(f"{i =:4d}, {t:.6g}/{tEnd:.6g}")
      if (i % outf) == 0:
          p.write_mig(outputdir, t)
      forces = np.pi * p.r() ** 2 * rho * g
      time_integration._advance_particles(p, forces, dt, 1, 1e-6 * r)
      t += dt
      i += 1
  

Plot
----
.. code-block:: shell

 python3 -m migflow.plot.migplot output --actors particles

