... | ... | @@ -4,7 +4,7 @@ Photoelatic image: inverse method analysis |
|
|
## 1. Overview
|
|
|
The inverse method analysis provides a quantitative estimation of the force network inside a granular media made in a photoelastic material. With such analysis, the contact force's magnitude and orientation can be determined under some assumptions. The main idea of the inverse method is to generate a numerical photoelastic picture that matches the experimental one like illustrated below (the experimental picture comes from [the Matlab implementation of the method](https://github.com/jekollmer/PEGS)).
|
|
|
|
|
|
<img src="uploads/7e9d15a3fc844774cb987ea3fb530f0e/Force_measurement.svg" width="800">
|
|
|
<img src="/uploads/7e9d15a3fc844774cb987ea3fb530f0e/Force_measurement.svg" width="800">
|
|
|
|
|
|
|
|
|
To achieve a good matching an optimization procedure is used, this procedure required a numerical computation of the photoelastic signal in the granular media. For cylindrical particles, such an analytical expression of this signal can be obtained using the theory of elasticity and the theory of photoelasticity.
|
... | ... | @@ -31,8 +31,8 @@ The assumptions required to obtain the stress field inside the disk are: |
|
|
5. No body-forces (gravity...);
|
|
|
|
|
|
The stress field inside the disk is obtained by the superposition of:
|
|
|
1. A simple radial stress distribution for each load $`i`$: $`\sigma_{r_i} = -\frac{2f_{i}}{h\pi}\frac{\cos\theta_i}{r_i}`$;
|
|
|
2. An uniform tension to ensure a stress-free boundary: $`\sigma_{xx} = \sigma_{yy} = \sum_{i=1}^{N}\frac{f_{i}}{\pi h d}\cos\left(\theta_{1,i}+\theta_{2,i}\right)`$.
|
|
|
1. A simple radial stress distribution for each load \(i\): \(\sigma_{r_i} = -\frac{2f_{i}}{h\pi}\frac{\cos\theta_i}{r_i}\);
|
|
|
2. An uniform tension to ensure a stress-free boundary: \(\sigma_{xx} = \sigma_{yy} = \sum_{i=1}^{N}\frac{f_{i}}{\pi h d}\cos\left(\theta_{1,i}+\theta_{2,i}\right)\).
|
|
|
|
|
|
Of course for the summation, the stress tensors have to be expressed in the same coordinates system.
|
|
|
|
... | ... | @@ -42,9 +42,9 @@ In a [circular polariscope](/reflection-photoelasticity#21-photoelasticity-in-th |
|
|
|
|
|
From the position and dimensions of the disk, it is easy to construct a possible contact network based on geometrical considerations. A contact occurs if the distance between the two disk centers is lower than the sum of the two radii increased by a small tolerance to overcome imprecision from the detection step. This provides the location of all possible contact forces acting on each disk, but it includes some false contacts (due to the tolerance) or _non-force-bearing_ (which does not transmit any load).
|
|
|
|
|
|
A correct result of the optimization procedure required a good initial force distribution. The [gradient analysis](/gradient-analysis) provides the $`G^2`$ value which is proportional to the sum of the contact force magnitude $`\sum_i|F_i|`$ acting on the disk. To provide the initial guess of force distribution, two methods are available depending on the quality of the pictures.
|
|
|
A correct result of the optimization procedure required a good initial force distribution. The [gradient analysis](/gradient-analysis) provides the \(G^2\) value which is proportional to the sum of the contact force magnitude \(\sum_i|F_i|\) acting on the disk. To provide the initial guess of force distribution, two methods are available depending on the quality of the pictures.
|
|
|
|
|
|
__In the case of high-quality pictures__ (resolution and contrast), a $`G^2_i`$ at each contact $`i`$ can be computed using a small region of interest located near the contact point. The false and _non-force-bearing_ contacts are eliminated if their $`G^2_i`$ value is lower than a threshold value on both contacting disks. The sum of the contact force magnitude is then distributed on the remaining contacts in proportion to the value of $`G^2_i`$. [1]
|
|
|
__In the case of high-quality pictures__ (resolution and contrast), a \(G^2_i\) at each contact \(i\) can be computed using a small region of interest located near the contact point. The false and _non-force-bearing_ contacts are eliminated if their \(G^2_i\) value is lower than a threshold value on both contacting disks. The sum of the contact force magnitude is then distributed on the remaining contacts in proportion to the value of \(G^2_i\). [1]
|
|
|
|
|
|
__In the case of low-quality pictures__, each possible combination of active contact is provided to the optimization procedure. The sum of the contact force magnitude is equally distributed over the active contacts. This approach is slower than the previous one and only gives good results when the contact force are almost equally distributed over the active contacts. To overcome the second drawback the fitting procedure will consider all the granular media to propagate the fitted forces on the adjacent disk when the optimization procedure succeeds (PhD Thesis of O. Lantsoght, Université Catholique de Louvain, 2019).
|
|
|
|
... | ... | @@ -54,7 +54,7 @@ __In the case of low-quality pictures__, each possible combination of active con |
|
|
|
|
|
The optimization procedure uses the force's magnitude and orientation as parameters and tries to minimize the difference between the computed photoelastic signal and the picture of the disk. Different methods exist to compare two images, among them we have the Mean Squared Error (MSE) or the Structural Similarity Index (SSIM) [3,4]. The main advantage of the MSE is the speed of computation but it considers the image globally. On the opposite, the SSIM is slower but considers the region of the image to take into account its structure. Presently, for this application, both comparison criteria are similar, they are globally good but in some cases failed to detect that the images do not match at all.
|
|
|
|
|
|
> <img src="uploads/46001dc532cf60a70fe7e5614e066855/SSIM_MSE.svg" width="400">
|
|
|
> <img src="/uploads/46001dc532cf60a70fe7e5614e066855/SSIM_MSE.svg" width="400">
|
|
|
>
|
|
|
> Comparison of SSIM and MSE index value on some particles (from Fig.1). The SSIM value range is between 0 and 1. The optimization process is a minimization process so we compute 1-SSIM.
|
|
|
|
... | ... | @@ -68,26 +68,26 @@ To improve the initial guess of the forces provided to the optimizer, the forces |
|
|
|
|
|
> The forces must be identified in order to reproduce the photoelastic signal of Fig.2. Figure 3 shows the gradient value at each point of the particles, the G2 value on each particle, and the sum of the contact force magnitude.
|
|
|
>
|
|
|
> <img src="uploads/1bc4999fead9eee0109faeb2ed086c96/Propa_1_2.svg" width="800">
|
|
|
> <img src="/uploads/1bc4999fead9eee0109faeb2ed086c96/Propa_1_2.svg" width="800">
|
|
|
>
|
|
|
> The optimizer is called for all particles with all combinations of active contact. The sum of the contact force magnitude is evenly distributed on all active contact. Figure 4 give all possible combination for particle number 1 (top left), the photoelastic signal obtained after optimization, and the SSIM value. The Fig.5 shows the result of this first run of optimization, the solved particles are the green ones.
|
|
|
>
|
|
|
> <img src="uploads/8725ddd3b0c2adab6d6e51b2a9e719c0/algoForceJeuxHautGauche_TXT.svg" width="800">
|
|
|
> <img src="/uploads/8725ddd3b0c2adab6d6e51b2a9e719c0/algoForceJeuxHautGauche_TXT.svg" width="800">
|
|
|
>
|
|
|
> The next iteration of optimization will be done on particles number 5 and 7 as they both have a solved contact from the neighbors (particles 3 and 4). All possible combinations for particle number 5 are illustrated in Fig.6. note that the contact force value with particles 3 and 4 are the same for all combinations. After this step of iteration, the particle 5 is solved but not the 7 (Fig.7).
|
|
|
>
|
|
|
> <img src="uploads/c67491baa350db50f7252e0499c9b7b4/second_it.svg" width="800">
|
|
|
> <img src="/uploads/c67491baa350db50f7252e0499c9b7b4/second_it.svg" width="800">
|
|
|
>
|
|
|
> The next iteration of optimization will only be carried out on particles 6 and 7. The only unknowns are the contact forces with the walls. This step allows us to retrieve the photoelastic response of Fig.3.
|
|
|
>
|
|
|
> The following pictures show the same algorithm on a more complex example. The first picture is the experimental photoelastic response. The next pictures illustrate the iterations 1 to 5, after that the algorithm does not reach the specified SSIM value on any new particle and the optimization process ends.
|
|
|
>
|
|
|
> <img src="uploads/cadedc89a2659bbcb5d0d927d5767e85/Picture.png" width="250">
|
|
|
> <img src="uploads/093ba577f1853595e148c32d324fb0da/synthetic0.png" width="250">
|
|
|
> <img src="uploads/7ec3e96162cda9ad086b560426877142/synthetic1.png" width="250">
|
|
|
> <img src="uploads/5cf649d69e2a1180a8a59511b339ceb9/synthetic2.png" width="250">
|
|
|
> <img src="uploads/7e5f9fa31da2d739f34f372157221460/synthetic3.png" width="250">
|
|
|
> <img src="uploads/a0ea651b60494eb790cef8ef30a2dac9/synthetic4.png" width="250">
|
|
|
> <img src="/uploads/cadedc89a2659bbcb5d0d927d5767e85/Picture.png" width="250">
|
|
|
> <img src="/uploads/093ba577f1853595e148c32d324fb0da/synthetic0.png" width="250">
|
|
|
> <img src="/uploads/7ec3e96162cda9ad086b560426877142/synthetic1.png" width="250">
|
|
|
> <img src="/uploads/5cf649d69e2a1180a8a59511b339ceb9/synthetic2.png" width="250">
|
|
|
> <img src="/uploads/7e5f9fa31da2d739f34f372157221460/synthetic3.png" width="250">
|
|
|
> <img src="/uploads/a0ea651b60494eb790cef8ef30a2dac9/synthetic4.png" width="250">
|
|
|
|
|
|
|
|
|
## 3. Examples
|
... | ... | @@ -95,16 +95,16 @@ To improve the initial guess of the forces provided to the optimizer, the forces |
|
|
### Force propagation
|
|
|
The following pictures show the picture of an experiment (left) and the corresponding photoelastic signal (right). In this experiment, the frame is vertically vibrated while the blender (stadium black shape at the center) oscillates by rotating around its upper extremity (that does not move vertically). In this experiment a circular bright-field polariscope is used, this is why the background is white.
|
|
|
|
|
|
<img src="uploads/6cb8360723ac75a52b43600e86d74bde/Wiki_Global.png" width="250">
|
|
|
<img src="uploads/31712778bfe2de3263c13d30a0b615cc/PhotoElast.png" width="250">
|
|
|
<img src="/uploads/6cb8360723ac75a52b43600e86d74bde/Wiki_Global.png" width="250">
|
|
|
<img src="/uploads/31712778bfe2de3263c13d30a0b615cc/PhotoElast.png" width="250">
|
|
|
|
|
|
The quality of the picture only allows us to correctly identify good enough initial forces to converge the photoelastic signal on disks with simple force distribution (top left figure). The force propagation procedure allows to reach convergence on more disks with more complex loading.
|
|
|
|
|
|
<img src="uploads/d5fa3d72a01cb149cfafe889ed594540/synthetic0.png" width="250">
|
|
|
<img src="uploads/ffb7930bcf8f1ab5fc19bebfb92726ad/synthetic1.png" width="250">
|
|
|
<img src="/uploads/d5fa3d72a01cb149cfafe889ed594540/synthetic0.png" width="250">
|
|
|
<img src="/uploads/ffb7930bcf8f1ab5fc19bebfb92726ad/synthetic1.png" width="250">
|
|
|
|
|
|
<img src="uploads/381567fa79c5355828309dfb539bbcae/synthetic3.png" width="250">
|
|
|
<img src="uploads/1f5ec507b8fef180cf598736c1f02c26/synthetic6.png" width="250">
|
|
|
<img src="/uploads/381567fa79c5355828309dfb539bbcae/synthetic3.png" width="250">
|
|
|
<img src="/uploads/1f5ec507b8fef180cf598736c1f02c26/synthetic6.png" width="250">
|
|
|
|
|
|
Eventually, the procedure did not succeed to reach convergences on all disks (bottom right figure) for reasons related to picture quality as well as assumption violation:
|
|
|
* Some disks show so close fringes that the picture of the signal is just gray, leading to completely wrong forces estimation;
|
... | ... | @@ -131,4 +131,4 @@ Photo-elastic Grain Solver (PEGS): |
|
|
|
|
|
|
|
|
|
|
|
[<go back to home](/home) |
|
|
\ No newline at end of file |
|
|
[<go back to home](/) |