cubo

.gitignore

 *.aux
 *.log
 *.out
-problem.ppl
-M_*.*
-MB_*.*
-pipe-*.png
-sigma*.*
-ur.*
-*-vs-*.*
-distorted*.png
-pipe-linearized.html
-pipe-linearized.tex
-pipe-linearized.pdf
+*.msh
+*.vtk
 nafems4.html
 nafems4.tex
 nafems4.pdf

cube.fin

+MESH FILE_PATH cube.msh
+
+E = 200e3
+nu = 0.3
+
+PHYSICAL_ENTITY left   BC u=0
+PHYSICAL_ENTITY front  BC v=0
+PHYSICAL_ENTITY bottom BC w=0
+
+# PHYSICAL_ENTITY right  BC tx=-1
+# PHYSICAL_ENTITY back   BC ty=2
+# PHYSICAL_ENTITY top    BC tz=3
+# 
+# PHYSICAL_ENTITY right  BC ty=1 tz=1
+# PHYSICAL_ENTITY back   BC tx=1 tz=1
+# PHYSICAL_ENTITY top    BC tx=1 ty=1
+
+PHYSICAL_ENTITY right  BC tx=-1 ty=1 tz=1
+PHYSICAL_ENTITY back   BC ty=2  tx=1 tz=1
+PHYSICAL_ENTITY top    BC tz=3  tx=1 ty=1
+
+
+FINO_SOLVER PROGRESS_ASCII
+FINO_STEP
+
+MESH_POST FILE_PATH cube-out3.msh VECTOR u v w sigma sigma1 sigma2 sigma3
+

cube.geo

+//+
+SetFactory("OpenCASCADE");
+//+
+a = 1;
+n = 20;
+Box(1) = {-0.5*a, -0.5*a, -0.5*a, a, a, a};
+
+Mesh.Algorithm = 6;
+Mesh.ElementOrder = 1;
+
+Physical Surface("left") = {1};
+Physical Surface("right") = {2};
+Physical Surface("front") = {3};
+Physical Surface("back") = {4};
+Physical Surface("bottom") = {5};
+Physical Surface("top") = {6};
+Physical Volume("bulk") = {1};
+
+
+Transfinite Line "*" = n;
+Transfinite Surface "*";
+Transfinite Volume "*";
+
+Mesh.RecombineAll = 1;
+Mesh.Recombine3DAll = 1;
+Mesh.Recombine3DLevel = 2;
+
+Mesh.Optimize = 1;
+

nafems4.md

 
 Just for the sake of it, let us jump out of our nuclear piping problem and step back into the general finite-element theory ground (remember we were going to jump back and forth). We have a linear elastic problem (well, we still do not have it but we will in\ [@sec:break]), so we might exploit the problem’s linearity and compute all the effects separately and them sum them up to obtain the whole combination, right? Well, not so much. Just like most of the time when we want to weight two masses we can sum their weights individually to obtain the same value we might get when putting both of them into a scale, most of the time we might just sum partial finite-element results. If we individually weighted two protons and two neutrons on the one hand and one $\alpha$ particle (i.e. a ⁴He nucleus) we would definitely not get the same result.
 
+divert(-1)
 Let us both (i.e. you and me) make an experiment. Grab a FEM program of your choice, get a square-sectioned beam of any size and length, fix one of the ends and put an uniform vertical load in the top surface. In my case, I have a square section of 1mm\ $\times$\ 1mm and a cantilever length of 10mm. A vertical load of 10N uniformly distributed in the top surface gives a maximum displacement of  and a maximum Von Mises stress of  
 
 https://caeplex.com/p?9c1
 podemos sumar los desplazamientos y los tensores de tensiones, pero ojo que las tensiones principales no son lineales a la suma (si a la multiplicacion)
 
 
+octave:32> A = rand(3); A = A*A'
+A =
+
+   2.08711   1.40929   1.31108
+   1.40929   1.32462   0.57570
+   1.31108   0.57570   1.09657
+
+octave:33> B = rand(3); B = B*B'
+B =
+
+   1.02619   0.73457   0.56903
+   0.73457   0.53386   0.37772
+   0.56903   0.37772   0.53141
+
+octave:34> eig(A)+eig(B)
+ans =
+
+   0.0075113
+   0.8248395
+   5.7674016
+
+octave:35> eig(A+B)
+ans =
+
+   0.049508
+   0.782990
+   5.767255
+
+octave:36> 
+
+divert(0)
+
+
+
 
 cantilever beam, principal stresses, linearity of von mises
 

principal.was

+sx1 = 1
+sy1 = 2
+sz1 = 3
+txy1 = 4
+tyz1 = 5
+tzx1 = 6
+
+sx2 = 7
+sy2 = 8
+sz2 = 9
+txy2 = 10
+tyz2 = 11
+tzx2 = 12
+
+IF $1=1
+sigmax = sx1
+sigmay = sy1
+sigmaz = sz1
+tauxy = txy1
+tauyz = tyz1
+tauzx = tzx1
+ENDIF
+IF $1=2
+sigmax = sx2
+sigmay = sy2
+sigmaz = sz2
+tauxy = txy2
+tauyz = tyz2
+tauzx = tzx2
+ENDIF
+IF $1=3
+sigmax = sx1 +  sx2
+sigmay = sy1 +  sy2
+sigmaz = sz1 +  sz2
+tauxy = txy1 + txy2
+tauyz = tyz1 + tyz2
+tauzx = tzx1 + tzx2
+ENDIF
+
+I1 = sigmax + sigmay + sigmaz
+I2 = sigmax*sigmay + sigmay*sigmaz + sigmaz*sigmax - tauxy^2 - tauyz^2 - tauzx^2
+I3 = sigmax*sigmay*sigmaz - sigmax*tauyz^2 - sigmay*tauzx^2 - sigmaz*tauxy^2 + 2*tauxy*tauyz*tauzx
+
+c1 = sqrt(abs(I1^2 - 3*I2))
+phi = 1.0/3.0 * acos((2.0*I1^3 - 9.0*I1*I2 + 27.0*I3)/(2.0*c1^3))
+c2 = I1/3.0
+c3 = 2.0/3.0 * c1
+
+sigma1 = c2 + c3 * cos(phi)
+sigma2 = c2 + c3 * cos(phi - 2.0*pi/3.0)
+sigma3 = c2 + c3 * cos(phi - 4.0*pi/3.0)
+
+PRINT I1 I2 I3 phi sigma1 sigma2 sigma3