diff --git a/examples/Freefem/validation/3D/3D_distri/beam3d_tet.msh b/examples/Freefem/validation/3D/3D_distri/beam3d_tet.msh
new file mode 100644
index 00000000..a6239eff
--- /dev/null
+++ b/examples/Freefem/validation/3D/3D_distri/beam3d_tet.msh
@@ -0,0 +1,588 @@
+$MeshFormat
+2.2 0 8
+$EndMeshFormat
+$PhysicalNames
+7
+2 1 "Fixed"
+2 2 "Bottom"
+2 3 "Right"
+2 4 "Top"
+2 5 "Front"
+2 6 "Back"
+3 7 "Beam"
+$EndPhysicalNames
+$Nodes
+105
+1 0 0 0
+2 1 0 0
+3 1 0.2 0
+4 0 0.2 0
+5 0 0 0.2
+6 1 0 0.2
+7 1 0.2 0.2
+8 0 0.2 0.2
+9 0.09999999999981413 0 0
+10 0.1999999999995569 0 0
+11 0.299999999999265 0 0
+12 0.3999999999989731 0 0
+13 0.4999999999986921 0 0
+14 0.599999999998945 0 0
+15 0.6999999999992088 0 0
+16 0.7999999999994725 0 0
+17 0.8999999999997362 0 0
+18 1 0.0999999999997371 0
+19 0.8999999999995836 0.2 0
+20 0.8 0.2 0
+21 0.7000000000006938 0.2 0
+22 0.6000000000013874 0.2 0
+23 0.5000000000020595 0.2 0
+24 0.400000000001665 0.2 0
+25 0.3000000000012487 0.2 0
+26 0.2000000000008325 0.2 0
+27 0.1000000000004162 0.2 0
+28 0 0.100000000000274 0
+29 0.09999999999981413 0 0.2
+30 0.1999999999995569 0 0.2
+31 0.299999999999265 0 0.2
+32 0.3999999999989731 0 0.2
+33 0.4999999999986921 0 0.2
+34 0.599999999998945 0 0.2
+35 0.6999999999992088 0 0.2
+36 0.7999999999994725 0 0.2
+37 0.8999999999997362 0 0.2
+38 1 0.0999999999997371 0.2
+39 0.8999999999995836 0.2 0.2
+40 0.8 0.2 0.2
+41 0.7000000000006938 0.2 0.2
+42 0.6000000000013874 0.2 0.2
+43 0.5000000000020595 0.2 0.2
+44 0.400000000001665 0.2 0.2
+45 0.3000000000012487 0.2 0.2
+46 0.2000000000008325 0.2 0.2
+47 0.1000000000004162 0.2 0.2
+48 0 0.100000000000274 0.2
+49 0 0 0.0999999999997371
+50 1 0 0.0999999999997371
+51 1 0.2 0.0999999999997371
+52 0 0.2 0.0999999999997371
+53 0.1000000000001152 0.1000000000002203 0
+54 0.2000000000001947 0.1000000000001666 0
+55 0.3000000000002568 0.1000000000001129 0
+56 0.4000000000003191 0.1000000000000592 0
+57 0.5000000000003757 0.1000000000000056 0
+58 0.6000000000001662 0.09999999999995186 0
+59 0.6999999999999513 0.09999999999989816 0
+60 0.7999999999997363 0.09999999999984448 0
+61 0.89999999999966 0.09999999999979081 0
+62 0.1000000000001152 0.1000000000002203 0.2
+63 0.2000000000001947 0.1000000000001666 0.2
+64 0.3000000000002568 0.1000000000001129 0.2
+65 0.4000000000003192 0.1000000000000592 0.2
+66 0.5000000000003757 0.1000000000000055 0.2
+67 0.6000000000001662 0.09999999999995185 0.2
+68 0.6999999999999513 0.09999999999989817 0.2
+69 0.7999999999997363 0.09999999999984449 0.2
+70 0.89999999999966 0.09999999999979081 0.2
+71 0.09999999999981413 0 0.0999999999997371
+72 0.1999999999995569 0 0.0999999999997371
+73 0.299999999999265 0 0.0999999999997371
+74 0.3999999999989731 0 0.0999999999997371
+75 0.4999999999986922 0 0.0999999999997371
+76 0.599999999998945 0 0.0999999999997371
+77 0.6999999999992088 0 0.0999999999997371
+78 0.7999999999994725 0 0.0999999999997371
+79 0.8999999999997362 0 0.0999999999997371
+80 0.8999999999995836 0.2 0.0999999999997371
+81 0.8 0.2 0.0999999999997371
+82 0.7000000000006938 0.2 0.0999999999997371
+83 0.6000000000013874 0.2 0.0999999999997371
+84 0.5000000000020595 0.2 0.0999999999997371
+85 0.400000000001665 0.2 0.0999999999997371
+86 0.3000000000012487 0.2 0.0999999999997371
+87 0.2000000000008325 0.2 0.0999999999997371
+88 0.1000000000004162 0.2 0.0999999999997371
+89 1 0.0999999999997371 0.0999999999997371
+90 0 0.100000000000274 0.0999999999997371
+91 0.2589961673760527 0.1674090681079447 0.02458352422436767
+92 0.4574915415816455 0.1498224237280481 0.0569962371169359
+93 0.3496440916187388 0.04605589271472998 0.04587953831718808
+94 0.1499999999999093 0.05117268292715365 0.1488273170728393
+95 0.8553903463552315 0.1617363547396499 0.161736354739664
+96 0.5499999999993299 0.05117268292752882 0.05117268292760564
+97 0.5439242897806678 0.1579084165060936 0.04687715009988654
+98 0.7647727272728269 0.1574476117266997 0.0483567026358275
+99 0.6589961673756821 0.1674090681076312 0.02458352422439182
+100 0.8499999999996345 0.05117268292775361 0.05117268292778051
+101 0.7625447544532288 0.05117268292776111 0.05117268292778123
+102 0.8314182914960002 0.1028377246782176 0.1020179369065382
+103 0.5499999999988187 0.1255389831945015 0.1255389831954118
+104 0.3545791141791101 0.1071281928871111 0.1048763752863662
+105 0.1500000000001687 0.125538983192814 0.07446101680704631
+$EndNodes
+$Elements
+464
+1 2 2 5 1 53 9 1
+2 2 2 5 1 28 53 1
+3 2 2 5 1 18 2 17
+4 2 2 5 1 3 18 61
+5 2 2 5 1 19 3 61
+6 2 2 5 1 4 27 28
+7 2 2 5 1 54 10 9
+8 2 2 5 1 53 54 9
+9 2 2 5 1 55 11 10
+10 2 2 5 1 54 55 10
+11 2 2 5 1 56 12 11
+12 2 2 5 1 55 56 11
+13 2 2 5 1 57 13 12
+14 2 2 5 1 56 57 12
+15 2 2 5 1 58 14 13
+16 2 2 5 1 57 58 13
+17 2 2 5 1 59 15 14
+18 2 2 5 1 58 59 14
+19 2 2 5 1 60 16 15
+20 2 2 5 1 59 60 15
+21 2 2 5 1 61 17 16
+22 2 2 5 1 60 61 16
+23 2 2 5 1 61 18 17
+24 2 2 5 1 20 19 60
+25 2 2 5 1 19 61 60
+26 2 2 5 1 21 20 59
+27 2 2 5 1 20 60 59
+28 2 2 5 1 22 21 58
+29 2 2 5 1 21 59 58
+30 2 2 5 1 23 22 57
+31 2 2 5 1 22 58 57
+32 2 2 5 1 24 23 56
+33 2 2 5 1 23 57 56
+34 2 2 5 1 25 24 55
+35 2 2 5 1 24 56 55
+36 2 2 5 1 26 25 54
+37 2 2 5 1 25 55 54
+38 2 2 5 1 27 26 53
+39 2 2 5 1 26 54 53
+40 2 2 5 1 27 53 28
+41 2 2 6 2 5 29 48
+42 2 2 6 2 37 6 70
+43 2 2 6 2 70 6 38
+44 2 2 6 2 39 38 7
+45 2 2 6 2 8 62 47
+46 2 2 6 2 48 62 8
+47 2 2 6 2 29 30 62
+48 2 2 6 2 48 29 62
+49 2 2 6 2 30 31 63
+50 2 2 6 2 62 30 63
+51 2 2 6 2 31 32 64
+52 2 2 6 2 63 31 64
+53 2 2 6 2 32 33 65
+54 2 2 6 2 64 32 65
+55 2 2 6 2 33 34 66
+56 2 2 6 2 65 33 66
+57 2 2 6 2 34 35 67
+58 2 2 6 2 66 34 67
+59 2 2 6 2 35 36 68
+60 2 2 6 2 67 35 68
+61 2 2 6 2 36 37 69
+62 2 2 6 2 68 36 69
+63 2 2 6 2 69 37 70
+64 2 2 6 2 70 38 39
+65 2 2 6 2 40 70 39
+66 2 2 6 2 41 69 40
+67 2 2 6 2 69 70 40
+68 2 2 6 2 42 68 41
+69 2 2 6 2 68 69 41
+70 2 2 6 2 43 67 42
+71 2 2 6 2 67 68 42
+72 2 2 6 2 44 66 43
+73 2 2 6 2 66 67 43
+74 2 2 6 2 45 65 44
+75 2 2 6 2 65 66 44
+76 2 2 6 2 46 64 45
+77 2 2 6 2 64 65 45
+78 2 2 6 2 47 63 46
+79 2 2 6 2 63 64 46
+80 2 2 6 2 62 63 47
+81 2 2 2 3 1 9 49
+82 2 2 2 3 17 2 79
+83 2 2 2 3 79 2 50
+84 2 2 2 3 5 71 29
+85 2 2 2 3 49 71 5
+86 2 2 2 3 37 50 6
+87 2 2 2 3 9 10 71
+88 2 2 2 3 49 9 71
+89 2 2 2 3 10 11 72
+90 2 2 2 3 71 10 72
+91 2 2 2 3 11 12 73
+92 2 2 2 3 72 11 73
+93 2 2 2 3 12 13 74
+94 2 2 2 3 73 12 74
+95 2 2 2 3 13 14 75
+96 2 2 2 3 74 13 75
+97 2 2 2 3 14 15 76
+98 2 2 2 3 75 14 76
+99 2 2 2 3 15 16 77
+100 2 2 2 3 76 15 77
+101 2 2 2 3 16 17 78
+102 2 2 2 3 77 16 78
+103 2 2 2 3 78 17 79
+104 2 2 2 3 29 72 30
+105 2 2 2 3 71 72 29
+106 2 2 2 3 30 73 31
+107 2 2 2 3 72 73 30
+108 2 2 2 3 31 74 32
+109 2 2 2 3 73 74 31
+110 2 2 2 3 32 75 33
+111 2 2 2 3 74 75 32
+112 2 2 2 3 33 76 34
+113 2 2 2 3 75 76 33
+114 2 2 2 3 34 77 35
+115 2 2 2 3 76 77 34
+116 2 2 2 3 35 78 36
+117 2 2 2 3 77 78 35
+118 2 2 2 3 36 79 37
+119 2 2 2 3 78 79 36
+120 2 2 2 3 79 50 37
+121 2 2 4 4 3 19 51
+122 2 2 4 4 27 4 88
+123 2 2 4 4 88 4 52
+124 2 2 4 4 7 80 39
+125 2 2 4 4 51 80 7
+126 2 2 4 4 47 52 8
+127 2 2 4 4 19 20 80
+128 2 2 4 4 51 19 80
+129 2 2 4 4 20 21 81
+130 2 2 4 4 80 20 81
+131 2 2 4 4 21 22 82
+132 2 2 4 4 81 21 82
+133 2 2 4 4 22 23 83
+134 2 2 4 4 82 22 83
+135 2 2 4 4 23 24 84
+136 2 2 4 4 83 23 84
+137 2 2 4 4 24 25 85
+138 2 2 4 4 84 24 85
+139 2 2 4 4 25 26 86
+140 2 2 4 4 85 25 86
+141 2 2 4 4 26 27 87
+142 2 2 4 4 86 26 87
+143 2 2 4 4 87 27 88
+144 2 2 4 4 39 81 40
+145 2 2 4 4 80 81 39
+146 2 2 4 4 40 82 41
+147 2 2 4 4 81 82 40
+148 2 2 4 4 41 83 42
+149 2 2 4 4 82 83 41
+150 2 2 4 4 42 84 43
+151 2 2 4 4 83 84 42
+152 2 2 4 4 43 85 44
+153 2 2 4 4 84 85 43
+154 2 2 4 4 44 86 45
+155 2 2 4 4 85 86 44
+156 2 2 4 4 45 87 46
+157 2 2 4 4 86 87 45
+158 2 2 4 4 46 88 47
+159 2 2 4 4 87 88 46
+160 2 2 4 4 88 52 47
+161 2 2 3 5 2 18 50
+162 2 2 3 5 18 3 89
+163 2 2 3 5 89 3 51
+164 2 2 3 5 6 89 38
+165 2 2 3 5 50 89 6
+166 2 2 3 5 38 51 7
+167 2 2 3 5 50 18 89
+168 2 2 3 5 89 51 38
+169 2 2 1 6 1 49 28
+170 2 2 1 6 28 90 4
+171 2 2 1 6 4 90 52
+172 2 2 1 6 90 5 48
+173 2 2 1 6 49 5 90
+174 2 2 1 6 52 48 8
+175 2 2 1 6 28 49 90
+176 2 2 1 6 90 48 52
+177 4 2 7 1 59 99 58 103
+178 4 2 7 1 60 98 59 102
+179 4 2 7 1 55 104 91 105
+180 4 2 7 1 55 91 54 105
+181 4 2 7 1 62 30 29 94
+182 4 2 7 1 63 30 62 94
+183 4 2 7 1 84 23 83 97
+184 4 2 7 1 13 58 57 96
+185 4 2 7 1 40 81 39 95
+186 4 2 7 1 22 83 23 97
+187 4 2 7 1 81 80 39 95
+188 4 2 7 1 13 14 58 96
+189 4 2 7 1 55 11 56 93
+190 4 2 7 1 23 56 57 92
+191 4 2 7 1 12 56 11 93
+192 4 2 7 1 24 56 23 92
+193 4 2 7 1 16 17 61 100
+194 4 2 7 1 60 16 61 100
+195 4 2 7 1 55 25 54 91
+196 4 2 7 1 25 26 54 91
+197 4 2 7 1 59 60 20 98
+198 4 2 7 1 21 58 59 99
+199 4 2 7 1 77 16 15 101
+200 4 2 7 1 78 16 77 101
+201 4 2 7 1 21 59 20 98
+202 4 2 7 1 22 58 21 99
+203 4 2 7 1 59 101 60 102
+204 4 2 7 1 78 101 77 102
+205 4 2 7 1 91 104 87 105
+206 4 2 7 1 92 103 85 104
+207 4 2 7 1 82 98 81 102
+208 4 2 7 1 86 87 91 104
+209 4 2 7 1 85 92 84 103
+210 4 2 7 1 59 58 76 103
+211 4 2 7 1 74 57 56 104
+212 4 2 7 1 57 75 103 104
+213 4 2 7 1 57 74 75 104
+214 4 2 7 1 55 73 104 105
+215 4 2 7 1 68 77 67 103
+216 4 2 7 1 68 78 77 102
+217 4 2 7 1 76 67 77 103
+218 4 2 7 1 69 78 68 102
+219 4 2 7 1 74 65 75 104
+220 4 2 7 1 63 64 73 104
+221 4 2 7 1 55 54 72 105
+222 4 2 7 1 55 72 73 105
+223 4 2 7 1 87 104 63 105
+224 4 2 7 1 83 68 67 103
+225 4 2 7 1 82 69 68 102
+226 4 2 7 1 69 82 81 102
+227 4 2 7 1 63 87 64 104
+228 4 2 7 1 86 64 87 104
+229 4 2 7 1 84 66 85 103
+230 4 2 7 1 85 24 84 92
+231 4 2 7 1 86 87 26 91
+232 4 2 7 1 26 25 86 91
+233 4 2 7 1 84 24 23 92
+234 4 2 7 1 83 22 82 99
+235 4 2 7 1 82 21 81 98
+236 4 2 7 1 20 81 21 98
+237 4 2 7 1 69 40 70 95
+238 4 2 7 1 21 82 22 99
+239 4 2 7 1 70 40 39 95
+240 4 2 7 1 22 57 58 97
+241 4 2 7 1 16 78 17 100
+242 4 2 7 1 79 17 78 100
+243 4 2 7 1 15 60 59 101
+244 4 2 7 1 15 16 60 101
+245 4 2 7 1 14 75 76 96
+246 4 2 7 1 74 12 73 93
+247 4 2 7 1 23 57 22 97
+248 4 2 7 1 13 75 14 96
+249 4 2 7 1 73 12 11 93
+250 4 2 7 1 30 72 29 94
+251 4 2 7 1 71 29 72 94
+252 4 2 7 1 60 101 78 102
+253 4 2 7 1 78 100 60 102
+254 4 2 7 1 75 96 57 103
+255 4 2 7 1 73 93 55 104
+256 4 2 7 1 63 72 94 105
+257 4 2 7 1 89 61 79 102
+258 4 2 7 1 71 90 62 105
+259 4 2 7 1 56 92 85 104
+260 4 2 7 1 83 58 99 103
+261 4 2 7 1 80 38 70 89
+262 4 2 7 1 38 80 51 89
+263 4 2 7 1 38 80 70 39
+264 4 2 7 1 62 90 71 29
+265 4 2 7 1 48 29 90 62
+266 4 2 7 1 74 57 75 13
+267 4 2 7 1 41 69 68 82
+268 4 2 7 1 70 6 89 38
+269 4 2 7 1 80 38 51 7
+270 4 2 7 1 6 70 89 37
+271 4 2 7 1 19 89 3 51
+272 4 2 7 1 62 8 52 48
+273 4 2 7 1 3 61 89 19
+274 4 2 7 1 31 73 63 64
+275 4 2 7 1 38 80 39 7
+276 4 2 7 1 5 71 49 90
+277 4 2 7 1 89 19 80 51
+278 4 2 7 1 19 89 80 61
+279 4 2 7 1 69 41 40 82
+280 4 2 7 1 5 71 90 29
+281 4 2 7 1 48 29 5 90
+282 4 2 7 1 61 3 89 18
+283 4 2 7 1 68 67 42 83
+284 4 2 7 1 37 89 6 50
+285 4 2 7 1 89 37 70 79
+286 4 2 7 1 82 40 69 81
+287 4 2 7 1 47 8 52 62
+288 4 2 7 1 62 52 90 48
+289 4 2 7 1 53 1 49 28
+290 4 2 7 1 30 73 63 31
+291 4 2 7 1 89 37 79 50
+292 4 2 7 1 73 30 63 72
+293 4 2 7 1 72 11 55 73
+294 4 2 7 1 49 53 90 71
+295 4 2 7 1 49 53 71 9
+296 4 2 7 1 83 41 68 82
+297 4 2 7 1 15 76 59 14
+298 4 2 7 1 53 28 49 90
+299 4 2 7 1 41 68 42 83
+300 4 2 7 1 57 74 12 13
+301 4 2 7 1 69 36 68 78
+302 4 2 7 1 1 53 49 9
+303 4 2 7 1 79 2 17 18
+304 4 2 7 1 52 90 88 62
+305 4 2 7 1 57 74 56 12
+306 4 2 7 1 63 87 46 64
+307 4 2 7 1 35 67 68 77
+308 4 2 7 1 27 28 88 4
+309 4 2 7 1 47 52 88 62
+310 4 2 7 1 10 72 55 54
+311 4 2 7 1 33 75 65 66
+312 4 2 7 1 79 18 89 50
+313 4 2 7 1 79 2 18 50
+314 4 2 7 1 11 72 55 10
+315 4 2 7 1 61 18 89 79
+316 4 2 7 1 79 61 18 17
+317 4 2 7 1 28 27 88 53
+318 4 2 7 1 85 44 65 66
+319 4 2 7 1 28 88 4 90
+320 4 2 7 1 35 78 68 36
+321 4 2 7 1 35 78 77 68
+322 4 2 7 1 14 76 59 58
+323 4 2 7 1 28 88 90 53
+324 4 2 7 1 90 88 4 52
+325 4 2 7 1 75 32 65 74
+326 4 2 7 1 45 87 86 64
+327 4 2 7 1 45 46 87 64
+328 4 2 7 1 77 34 67 76
+329 4 2 7 1 34 77 67 35
+330 4 2 7 1 43 66 85 84
+331 4 2 7 1 32 75 65 33
+332 4 2 7 1 85 44 66 43
+333 4 2 7 1 61 100 79 102
+334 4 2 7 1 58 96 76 103
+335 4 2 7 1 56 93 74 104
+336 4 2 7 1 62 94 71 105
+337 4 2 7 1 89 79 70 102
+338 4 2 7 1 89 80 61 102
+339 4 2 7 1 88 62 90 105
+340 4 2 7 1 71 53 90 105
+341 4 2 7 1 20 60 19 102
+342 4 2 7 1 19 60 61 102
+343 4 2 7 1 79 78 36 102
+344 4 2 7 1 36 37 79 102
+345 4 2 7 1 56 24 55 104
+346 4 2 7 1 55 24 25 104
+347 4 2 7 1 33 76 75 103
+348 4 2 7 1 74 73 31 104
+349 4 2 7 1 31 32 74 104
+350 4 2 7 1 33 34 76 103
+351 4 2 7 1 62 47 63 105
+352 4 2 7 1 46 63 47 105
+353 4 2 7 1 71 72 10 105
+354 4 2 7 1 9 71 10 105
+355 4 2 7 1 30 63 72 94
+356 4 2 7 1 78 60 16 101
+357 4 2 7 1 80 70 39 95
+358 4 2 7 1 11 55 73 93
+359 4 2 7 1 57 75 13 96
+360 4 2 7 1 60 78 16 100
+361 4 2 7 1 80 81 20 102
+362 4 2 7 1 69 70 37 102
+363 4 2 7 1 36 69 37 102
+364 4 2 7 1 20 19 80 102
+365 4 2 7 1 85 86 25 104
+366 4 2 7 1 33 66 34 103
+367 4 2 7 1 85 25 24 104
+368 4 2 7 1 31 64 32 104
+369 4 2 7 1 66 67 34 103
+370 4 2 7 1 32 64 65 104
+371 4 2 7 1 88 87 46 105
+372 4 2 7 1 46 47 88 105
+373 4 2 7 1 54 53 9 105
+374 4 2 7 1 9 10 54 105
+375 4 2 7 1 80 89 70 102
+376 4 2 7 1 88 90 53 105
+377 4 2 7 1 78 79 100 102
+378 4 2 7 1 75 76 96 103
+379 4 2 7 1 74 93 73 104
+380 4 2 7 1 71 94 72 105
+381 4 2 7 1 60 100 61 102
+382 4 2 7 1 57 96 58 103
+383 4 2 7 1 55 93 56 104
+384 4 2 7 1 62 63 94 105
+385 4 2 7 1 80 19 61 102
+386 4 2 7 1 36 78 69 102
+387 4 2 7 1 70 79 37 102
+388 4 2 7 1 56 85 24 104
+389 4 2 7 1 33 75 66 103
+390 4 2 7 1 73 64 31 104
+391 4 2 7 1 32 65 74 104
+392 4 2 7 1 67 76 34 103
+393 4 2 7 1 46 87 63 105
+394 4 2 7 1 9 53 71 105
+395 4 2 7 1 62 88 47 105
+396 4 2 7 1 10 72 54 105
+397 4 2 7 1 70 95 80 102
+398 4 2 7 1 83 97 58 103
+399 4 2 7 1 57 58 97 103
+400 4 2 7 1 71 62 29 94
+401 4 2 7 1 79 61 17 100
+402 4 2 7 1 22 83 58 99
+403 4 2 7 1 83 22 58 97
+404 4 2 7 1 74 56 12 93
+405 4 2 7 1 14 76 58 96
+406 4 2 7 1 85 56 24 92
+407 4 2 7 1 43 66 84 103
+408 4 2 7 1 40 69 81 102
+409 4 2 7 1 45 64 86 104
+410 4 2 7 1 42 83 67 103
+411 4 2 7 1 85 65 44 104
+412 4 2 7 1 53 27 88 105
+413 4 2 7 1 53 26 27 105
+414 4 2 7 1 53 54 26 105
+415 4 2 7 1 87 27 26 105
+416 4 2 7 1 87 88 27 105
+417 4 2 7 1 86 85 44 104
+418 4 2 7 1 44 45 86 104
+419 4 2 7 1 84 83 42 103
+420 4 2 7 1 42 43 84 103
+421 4 2 7 1 65 64 45 104
+422 4 2 7 1 44 65 45 104
+423 4 2 7 1 43 67 66 103
+424 4 2 7 1 42 67 43 103
+425 4 2 7 1 40 81 95 102
+426 4 2 7 1 87 26 91 105
+427 4 2 7 1 40 95 69 102
+428 4 2 7 1 54 91 26 105
+429 4 2 7 1 84 97 83 103
+430 4 2 7 1 69 95 70 102
+431 4 2 7 1 81 80 95 102
+432 4 2 7 1 105 73 63 72
+433 4 2 7 1 105 63 73 104
+434 4 2 7 1 83 68 99 82
+435 4 2 7 1 83 99 68 103
+436 4 2 7 1 59 68 99 103
+437 4 2 7 1 75 66 104 65
+438 4 2 7 1 75 104 66 103
+439 4 2 7 1 85 104 66 65
+440 4 2 7 1 85 66 104 103
+441 4 2 7 1 98 102 20 81
+442 4 2 7 1 20 102 98 60
+443 4 2 7 1 92 103 97 84
+444 4 2 7 1 97 103 92 57
+445 4 2 7 1 97 23 92 84
+446 4 2 7 1 92 23 97 57
+447 4 2 7 1 91 104 25 86
+448 4 2 7 1 25 104 91 55
+449 4 2 7 1 104 57 92 103
+450 4 2 7 1 104 92 57 56
+451 4 2 7 1 102 98 68 82
+452 4 2 7 1 102 68 98 59
+453 4 2 7 1 68 98 99 82
+454 4 2 7 1 68 99 98 59
+455 4 2 7 1 99 98 21 82
+456 4 2 7 1 99 21 98 59
+457 4 2 7 1 102 101 68 59
+458 4 2 7 1 102 68 101 77
+459 4 2 7 1 68 101 103 59
+460 4 2 7 1 68 103 101 77
+461 4 2 7 1 15 77 59 76
+462 4 2 7 1 15 59 77 101
+463 4 2 7 1 103 59 77 76
+464 4 2 7 1 103 77 59 101
+$EndElements
diff --git a/examples/Freefem/validation/3D/3D_distri/comparaison_script3d.py b/examples/Freefem/validation/3D/3D_distri/comparaison_script3d.py
new file mode 100644
index 00000000..c8eb6c2c
--- /dev/null
+++ b/examples/Freefem/validation/3D/3D_distri/comparaison_script3d.py
@@ -0,0 +1,130 @@
+"""
+3D Beam - Distributed Load on Top Face - Comparison File
+
+No analytical solution available for this case -> pure cross-validation
+SOFA vs FreeFEM, RMS computed separately on ux, uy, uz.
+"""
+import json
+import os
+import sys
+import numpy as np
+import matplotlib.pyplot as plt
+
+from sofa_beam3d_distributed import sofaRun
+from pyfreefem import FreeFemRunner
+
+
+def _rms(a, b):
+    """Discrete RMS, normalized by sqrt(n)."""
+    return np.linalg.norm(a - b) / np.sqrt(a.size)
+
+
+def _default_params_path():
+    return os.path.join(os.path.dirname(os.path.abspath(__file__)), "params_beam3d_distributed.json")
+
+
+def _default_mesh_path():
+    return os.path.join(os.path.dirname(os.path.abspath(__file__)), "beam3d_tet.msh")
+
+
+def _pair_by_coordinates(x_a, y_a, z_a, x_b, y_b, z_b, tol=1e-6): 
+    order_a = np.lexsort((z_a, y_a, x_a))
+    order_b = np.lexsort((z_b, y_b, x_b))
+    if not (np.allclose(x_a[order_a], x_b[order_b], atol=tol)
+            and np.allclose(y_a[order_a], y_b[order_b], atol=tol)
+            and np.allclose(z_a[order_a], z_b[order_b], atol=tol)):
+        raise ValueError("Node coordinates don't match between SOFA and FreeFEM meshes.")
+    perm = np.empty_like(order_b)
+    perm[order_b] = order_a
+    return perm
+
+
+if __name__ == "__main__":
+
+    config_file = sys.argv[1] if len(sys.argv) > 1 else _default_params_path()
+    with open(config_file) as f:
+        cfg = json.load(f)
+
+    q             = float(cfg["q"])
+    young_modulus = float(cfg["youngModulus"])
+    poisson_ratio = float(cfg["poissonRatio"])
+    mesh_file     = _default_mesh_path()
+
+    # Run FreeFEM 
+    runner = FreeFemRunner("freefem_beam3d_distributed.edp")
+    exports = runner.execute({
+        'q': q,
+        'youngModulus': young_modulus,
+        'poissonRatio': poisson_ratio,
+        'meshFile': mesh_file,
+    })
+    x_ff  = exports['xcoords']
+    y_ff  = exports['ycoords']
+    z_ff  = exports['zcoords']
+    ux_ff = exports['ux[]']
+    uy_ff = exports['uy[]']
+    uz_ff = exports['uz[]']
+
+    
+    pos0_sofa, u_sofa = sofaRun(mesh_file=mesh_file, q=q,
+                                 young_modulus=young_modulus,
+                                 poisson_ratio=poisson_ratio)
+    x_sofa, y_sofa, z_sofa = pos0_sofa[:, 0], pos0_sofa[:, 1], pos0_sofa[:, 2]
+    ux_sofa, uy_sofa, uz_sofa = u_sofa[:, 0], u_sofa[:, 1], u_sofa[:, 2]
+
+    
+    perm = _pair_by_coordinates(x_sofa, y_sofa, z_sofa, x_ff, y_ff, z_ff)
+    ux_ff_p = ux_ff[perm]
+    uy_ff_p = uy_ff[perm]
+    uz_ff_p = uz_ff[perm]
+
+    rms_ux = _rms(ux_sofa, ux_ff_p)
+    rms_uy = _rms(uy_sofa, uy_ff_p)
+    rms_uz = _rms(uz_sofa, uz_ff_p)
+
+    # --- Write results ---
+    os.makedirs("results", exist_ok=True)
+    with open("results/comparison_beam3d_distributed_results.txt", 'w') as f:
+        header = (f"{'x':>10}  {'y':>10}  {'z':>10}  "
+                  f"{'ux_sofa':>12}  {'ux_ff':>12}  "
+                  f"{'uy_sofa':>12}  {'uy_ff':>12}  "
+                  f"{'uz_sofa':>12}  {'uz_ff':>12}")
+        f.write(header + "\n")
+        f.write("-" * len(header) + "\n")
+        for x, y, z, uxs, uxf, uys, uyf, uzs, uzf in zip(
+                x_sofa, y_sofa, z_sofa, ux_sofa, ux_ff_p, uy_sofa, uy_ff_p, uz_sofa, uz_ff_p):
+            f.write(f"{x:10.4f}  {y:10.4f}  {z:10.4f}  "
+                    f"{uxs:12.6e}  {uxf:12.6e}  "
+                    f"{uys:12.6e}  {uyf:12.6e}  "
+                    f"{uzs:12.6e}  {uzf:12.6e}\n")
+
+        f.write("\n")
+        f.write("RMS norms (SOFA vs FreeFEM, discrete, sqrt(n)-normalized)\n")
+        f.write("-" * 55 + "\n")
+        f.write(f"  RMS_ux (SOFA vs FF) = {rms_ux:.6e}\n")
+        f.write(f"  RMS_uy (SOFA vs FF) = {rms_uy:.6e}\n")
+        f.write(f"  RMS_uz (SOFA vs FF) = {rms_uz:.6e}\n")
+
+    print(f"RMS_ux (SOFA vs FF) = {rms_ux:.6e}")
+    print(f"RMS_uy (SOFA vs FF) = {rms_uy:.6e}")
+    print(f"RMS_uz (SOFA vs FF) = {rms_uz:.6e}")
+ 
+    fig, axes = plt.subplots(1, 3, figsize=(15, 5))
+    fig.suptitle("3D Beam — Distributed Load on Top — SOFA vs FreeFEM (parity)", fontsize=14)
+
+    def _parity(ax, a, b, label):
+        ax.scatter(a, b, s=8, alpha=0.6)
+        lims = [min(a.min(), b.min()), max(a.max(), b.max())]
+        ax.plot(lims, lims, 'r--', linewidth=1)
+        ax.set_xlabel(f"{label}_sofa")
+        ax.set_ylabel(f"{label}_ff")
+        ax.set_title(label)
+        ax.set_aspect('equal')
+
+    _parity(axes[0], ux_sofa, ux_ff_p, "ux")
+    _parity(axes[1], uy_sofa, uy_ff_p, "uy")
+    _parity(axes[2], uz_sofa, uz_ff_p, "uz")
+    plt.tight_layout()
+    fig.savefig("results/comparison_beam3d_distributed_fields.png", dpi=150)
+    plt.close(fig)
+ 
\ No newline at end of file
diff --git a/examples/Freefem/validation/3D/3D_distri/freefem_beam3d_distributed.edp b/examples/Freefem/validation/3D/3D_distri/freefem_beam3d_distributed.edp
new file mode 100644
index 00000000..10a7a400
--- /dev/null
+++ b/examples/Freefem/validation/3D/3D_distri/freefem_beam3d_distributed.edp
@@ -0,0 +1,50 @@
+IMPORT "io.edp"
+load "gmsh"
+load "msh3"
+
+DEFAULT (q, 1000.0)
+DEFAULT (youngModulus, 1000.0)
+DEFAULT (poissonRatio, 0.3)
+DEFAULT (meshFile, "beam3d_tet.msh")
+
+real Q  = $q;
+real E  = $youngModulus;
+real nu = $poissonRatio;
+
+mesh3 Th = gmshload3("$meshFile");
+
+fespace Vh(Th, P1);
+Vh ux, uy, uz, vx, vy, vz;
+
+real lambda = E*nu / ((1.+nu)*(1.-2.*nu));
+real mu     = E / (2.*(1.+nu));
+
+
+problem Elasticity([ux, uy, uz], [vx, vy, vz]) =
+    int3d(Th)(
+          lambda*(dx(ux)+dy(uy)+dz(uz))*(dx(vx)+dy(vy)+dz(vz))
+        + 2.*mu*(dx(ux)*dx(vx) + dy(uy)*dy(vy) + dz(uz)*dz(vz))
+        + mu*(dy(ux)+dx(uy))*(dy(vx)+dx(vy))
+        + mu*(dz(ux)+dx(uz))*(dz(vx)+dx(vz))
+        + mu*(dz(uy)+dy(uz))*(dz(vy)+dy(vz))
+    )
+    - int2d(Th, 4)( -Q*vy )
+    + on(1, ux=0, uy=0, uz=0);
+
+Elasticity;
+
+exportArray(ux[]);
+exportArray(uy[]);
+exportArray(uz[]);
+
+real[int] xcoords(Th.nv);
+real[int] ycoords(Th.nv);
+real[int] zcoords(Th.nv);
+for (int i = 0; i < Th.nv; i++) {
+    xcoords[i] = Th(i).x;
+    ycoords[i] = Th(i).y;
+    zcoords[i] = Th(i).z;
+}
+exportArray(xcoords);
+exportArray(ycoords);
+exportArray(zcoords);
diff --git a/examples/Freefem/validation/3D/3D_distri/params_beam3d_distributed.json b/examples/Freefem/validation/3D/3D_distri/params_beam3d_distributed.json
new file mode 100644
index 00000000..817dc50d
--- /dev/null
+++ b/examples/Freefem/validation/3D/3D_distri/params_beam3d_distributed.json
@@ -0,0 +1,7 @@
+{
+
+"q": 1000,
+"youngModulus": 1000,
+"poissonRatio": 0.3
+
+}
\ No newline at end of file
diff --git a/examples/Freefem/validation/3D/3D_distri/sofa_beam3d_distributed.py b/examples/Freefem/validation/3D/3D_distri/sofa_beam3d_distributed.py
new file mode 100644
index 00000000..e7880dfb
--- /dev/null
+++ b/examples/Freefem/validation/3D/3D_distri/sofa_beam3d_distributed.py
@@ -0,0 +1,162 @@
+"""
+3D Beam Simulation - Distributed Load on Top Face - P1 Tetrahedra
+Cross-validation SOFA vs FreeFEM
+"""
+import json
+import os
+import sys
+import numpy as np
+import Sofa
+import Sofa.Core
+import Sofa.Simulation
+
+RESULTS_DIR = "results"
+
+
+def consistent_traction_forces(nodes, top_faces, q): 
+    N = len(nodes)
+    F = np.zeros((N, 3))
+    for tri in top_faces:
+        pts = nodes[tri, :]
+        v1 = pts[1] - pts[0]
+        v2 = pts[2] - pts[0]
+        area = 0.5 * np.linalg.norm(np.cross(v1, v2))
+        for nid in tri:
+            F[nid, 1] += -q * area / 3.0
+    return F
+
+
+def _default_mesh_path():
+    return os.path.join(os.path.dirname(os.path.abspath(__file__)), "beam3d_tet.msh")
+
+
+def _default_params_path():
+    return os.path.join(os.path.dirname(os.path.abspath(__file__)), "params_beam3d_distributed.json")
+
+
+def create_scene_args(rootNode, mesh_file, q, young_modulus, poisson_ratio, tol=1e-6):
+    requiredPlugins = [
+        "Elasticity",
+        "Sofa.Component.Constraint.Projective",
+        "Sofa.Component.IO.Mesh",
+        "Sofa.Component.LinearSolver.Direct",
+        "Sofa.Component.MechanicalLoad",
+        "Sofa.Component.ODESolver.Backward",
+        "Sofa.Component.StateContainer",
+        "Sofa.Component.Topology.Container.Dynamic",
+        "Sofa.Component.Visual",
+        "Sofa.GL.Component.Rendering3D",
+    ]
+    rootNode.addObject('RequiredPlugin', pluginName=requiredPlugins)
+    rootNode.addObject('DefaultAnimationLoop')
+    rootNode.addObject('VisualStyle', displayFlags=["showBehaviorModels", "showForceFields"])
+
+    template = "Vec3d"
+
+    with rootNode.addChild('Beam') as Beam:
+        Beam.addObject('NewtonRaphsonSolver'
+                    , name="newtonSolver"
+                    , printLog=True
+                    , maxNbIterationsNewton=30
+                    , absoluteResidualStoppingThreshold=1e-12)
+        Beam.addObject('SparseLDLSolver'
+                    , name="linearSolver"
+                    , template="CompressedRowSparseMatrixd")
+        Beam.addObject('StaticSolver'
+                    , name="staticSolver"
+                    , newtonSolver="@newtonSolver"
+                    , linearSolver="@linearSolver")
+
+        loader = Beam.addObject('MeshGmshLoader', name="loader", filename=mesh_file)
+
+        nodes = np.array(loader.position.value)
+        tets = np.array(loader.tetrahedra.value)
+        tris = np.array(loader.triangles.value)
+        N = len(nodes)
+
+        
+        x_min = nodes[:, 0].min()
+        y_max = nodes[:, 1].max()
+        fixed_idx = np.where(np.isclose(nodes[:, 0], x_min, atol=tol))[0].tolist()
+        top_mask = np.all(np.isclose(nodes[tris, 1], y_max, atol=tol), axis=1)
+        top_faces = tris[top_mask]
+
+        F_nodal = consistent_traction_forces(nodes, top_faces, q)
+        forces_list = F_nodal.tolist()
+
+        dofs = Beam.addObject('MechanicalObject'
+                            , name="dofs"
+                            , template=template
+                            , position="@loader.position"
+                            , showObject=True
+                            , showObjectScale=0.01)
+
+        Beam.addObject('TetrahedronSetTopologyContainer'
+                    , name="topology"
+                    , src="@loader")
+        Beam.addObject('TetrahedronSetTopologyModifier')
+
+        Beam.addObject('LinearSmallStrainFEMForceField'
+                    , name="FEM"
+                    , template=template
+                    , youngModulus=young_modulus
+                    , poissonRatio=poisson_ratio
+                    , topology="@topology")
+
+        Beam.addObject('FixedProjectiveConstraint'
+                    , name="dirichlet"
+                    , indices=fixed_idx)
+
+        Beam.addObject('ConstantForceField'
+                    , name="TopTraction"
+                    , indices=list(range(N))
+                    , forces=forces_list)
+
+    return rootNode, dofs, nodes.copy()
+
+
+def createScene(rootNode):
+    with open(_default_params_path()) as f:
+        cfg = json.load(f)
+    create_scene_args(rootNode
+                    , mesh_file=_default_mesh_path()
+                    , q=float(cfg["q"])
+                    , young_modulus=float(cfg["youngModulus"])
+                    , poisson_ratio=float(cfg["poissonRatio"]))
+    return rootNode
+
+
+def sofaRun(mesh_file, q, young_modulus, poisson_ratio):
+    root = Sofa.Core.Node("root")
+    _, dofs, pos0 = create_scene_args(root
+                                    , mesh_file=mesh_file
+                                    , q=q
+                                    , young_modulus=young_modulus
+                                    , poisson_ratio=poisson_ratio)
+    Sofa.Simulation.init(root)
+    Sofa.Simulation.animate(root, root.dt.value)
+
+    pos_final = np.array(dofs.position.toList())
+    u = pos_final[:, :3] - pos0[:, :3]
+    x0 = pos0[:, :3]
+
+    os.makedirs(RESULTS_DIR, exist_ok=True)
+    out_path = os.path.join(RESULTS_DIR, "sofa_beam3d_distributed_results.txt")
+    with open(out_path, 'w') as f:
+        f.write(f"{'x0':>12}  {'y0':>12}  {'z0':>12}  {'ux':>12}  {'uy':>12}  {'uz':>12}\n")
+        f.write("-" * 78 + "\n")
+        for (xi, yi, zi), (uxi, uyi, uzi) in zip(x0, u):
+            f.write(f"{xi:12.6f}  {yi:12.6f}  {zi:12.6f}  {uxi:12.6f}  {uyi:12.6f}  {uzi:12.6f}\n")
+
+    return x0, u
+
+
+if __name__ == "__main__":
+    config_file = sys.argv[1] if len(sys.argv) > 1 else _default_params_path()
+    with open(config_file) as f:
+        cfg = json.load(f)
+
+    sofaRun(mesh_file=_default_mesh_path()
+            , q=float(cfg["q"])
+            , young_modulus=float(cfg["youngModulus"])
+            , poisson_ratio=float(cfg["poissonRatio"]))
\ No newline at end of file