changes after angus call

ams.tex

@@ -0,0 +1 @@
+\usepackage{amstext}

meta.yaml

@@ -1,5 +1,5 @@
 ---
-title: The NAFEMS Benchmark Challenge Volume\ 1
+title: The NAFEMS Case Studies Volumes\ 1
 author: Angus Ramsay
 fontsize: 10pt
 lang: en-GB

nafems4-emf.md

@@ -9,13 +9,18 @@
 
 # CS4: Fatigue Assessment in Nuclear Piping
 
-## From college theory to an actual engineering problem
+## From college theory to an actual engineering problem {#sec:intro}
 
 First of all, please take this text as a written chat between you an me, i.e. an average engineer that has already taken the journey from college to performing actual engineering works using finite element analysis and has something to say about it. Picture yourself in a coffee bar, talking and discussing concepts and ideas with me. Maybe needing to go to a blackboard (or notepad?). Even using a tablet to illustrate some three-dimensional results. But always as a chat between colleagues.
 
 Please also note that I am not a mechanical engineer, although I shared many undergraduate courses with some of them. I am a nuclear engineer with a strong background in mathematics and computer programming. I went to college between 2002 and 2008. Probably a lot of things have changed since then---at least that is what these “millennial” guys and girls seem to be boasting about---but chances are we all studied solid mechanics and heat transfer with a teacher using a piece of chalk on a blackboard while we as students wrote down notes with pencils on paper sheets. And there is really not much that one can do with pencil and paper regarding mechanical analysis. Any actual case worth the time of an engineer needs to be more complex than an ideal canonical case with analytical solution.
 
-Whether you are a student or a seasoned engineer with many years of experience, you might recall from first year physics courses the introduction of the [simple pendulum](https://en.wikipedia.org/wiki/Pendulum) as a case study\ ([@fig:simple]). You learned that the period does not depend on the hanging mass because the weight and the inertia exactly cancelled each other. Also, that [Galileo](https://en.wikipedia.org/wiki/Galileo_Galilei) said (and [Newton proved](https://www.seamplex.com/wasora/realbook/real-012-mechanics.html)) that for small oscillations the period does not even depend on the amplitude. Someone showed you why it worked this way: because if\ $\sin \theta \approx \theta$ then the motion equations converge to an [harmonic oscillator](https://en.wikipedia.org/wiki/Harmonic_oscillator). It might have been a difficult subject for you back in those days when you were learning physics and calculus at the same time. You might have later studied the [Lagrangian](https://en.wikipedia.org/wiki/Lagrangian_mechanics) and even the [Hamiltonian](https://en.wikipedia.org/wiki/Hamiltonian_mechanics) formulations, added a [parametric excitation](https://en.wikipedia.org/wiki/Parametric_oscillator) and analysed the [chaotic double pendulum](https://www.seamplex.com/wasora/realbook/real-017-double-pendulum.html). But it was probably after college, say when you took your first son to a swing on a windy day ([@fig:hamaca]), that you were faced with a real pendulum worth your full attention. Yes, this is my personal story but it could have easily been yours as well. My point is that the very same distance between what I imagined studying a pendulum was and what I saw that day at the swing (namely that the period _does_ depend on the hanging mass) is the same distance between the mechanical problems studied in college and the actual cases encountered during a professional engineer’s lifetime ([@fig:pipes]). In this regard, I am referring only to technical issues. The part of dealing with clients, colleagues, bosses, etc. which is definitely not taught in engineering schools (you can get a heads up in business schools, but again it would be a theoretical pendulum) is way beyond the scope of both this article and my own capacities.
+Whether you are a student or a seasoned engineer with many years of experience, you might recall from first year physics courses the introduction of the [simple pendulum](https://en.wikipedia.org/wiki/Pendulum) as a case study\ ([@fig:simple]). You learned that the period does not depend on the hanging mass because the weight and the inertia exactly cancelled each other. Also, that [Galileo](https://en.wikipedia.org/wiki/Galileo_Galilei) said (and [Newton proved](https://seamplex.com/wasora/doc/realbook/012-mechanics/)) that for small oscillations the period does not even depend on the amplitude. Someone showed you why it worked this way: because if\ $\sin \theta \approx \theta$ then the motion equations converge to an [harmonic oscillator](https://en.wikipedia.org/wiki/Harmonic_oscillator). It might have been a difficult subject for you back in those days when you were learning physics and calculus at the same time.
+divert(-1)
+You might have later studied the [Lagrangian](https://en.wikipedia.org/wiki/Lagrangian_mechanics) and even the [Hamiltonian](https://en.wikipedia.org/wiki/Hamiltonian_mechanics) formulations, added a [parametric excitation](https://en.wikipedia.org/wiki/Parametric_oscillator) and analysed the [chaotic double pendulum](https://seamplex.com/wasora/doc/realbook/017-double-pendulum). divert(0)
+But it was probably after college, say when you took your first son to a swing on a windy day ([@fig:hamaca]), that you were faced with a real pendulum worth your full attention.
+dnl Yes, this is my personal story but it could have easily been yours as well.
+The very same distance between what I imagined studying a pendulum was and what I saw that day at the swing (namely that the period _does_ depend on the hanging mass) is the same distance between the mechanical problems studied in college and the actual cases encountered during a professional engineer’s lifetime ([@fig:pipes]). In this regard, I am referring only to technical issues. The part of dealing with clients, colleagues, bosses, etc. which is definitely not taught in engineering schools (you can get a heads up in business schools, but again it would be a theoretical pendulum) is way beyond the scope of both this article and my own capacities.
 
 ::::: {#fig:pendulum}
 ![Simple pendulum](simple.emf){#fig:simple width=35%}\  
@@ -28,17 +33,18 @@ A [simple pendulum from college physics courses](https://www.seamplex.com/wasora
 ![College pipe](infinite-pipe.emf){#fig:infinite-pipe width=40%}
 ![Real-life pipe](isometric.emf){#fig:isometric width=58%}
 
-An infinitely-long pressurised thick pipe as taught in college and an isometric drawing of a section of a real-life piping system
+An infinitely-long pressurised thick pipe as taught in college and an isometric drawing of a section of a real-life piping system. The details of the isometric are not important individually but they are included to emphasize the fact that a real engineering drawing is far more complex and has far more information than what a professor can draw in a class blackboard.
 :::::
 
+divert(-1)
 Like the pendulums above, we will be swinging back and forth between a case study about fatigue assessment of piping systems in a nuclear power plant and more generic topics related to finite elements and computational mechanics. These latter regressions will not remain just as abstract theoretical ideas. Not only will they be directly applicable to the development of the main case, but they will also apply to a great deal of other engineering problems tackled with the finite element method (and its cousins, [@sec:formulations]).
 
 \medskip
 
 Finite elements are like magic to me. I mean, I can follow the whole derivation of the equations, from the strong, weak and variational formulations of the equilibrium equations for the mechanical problem (or the energy conservation for heat transfer) down to the [algebraic multigrid](https://en.wikipedia.org/wiki/Multigrid_method) [preconditioner](https://en.wikipedia.org/wiki/Preconditioner) for the inversion of the [stiffness matrix](https://en.wikipedia.org/wiki/Stiffness_matrix) passing through [Sobolev spaces](https://en.wikipedia.org/wiki/Sobolev_space) and the [grid generation](https://en.wikipedia.org/wiki/Mesh_generation). Then I can sit down and program all these steps into a computer, including the [shape functions](https://en.wikipedia.org/wiki/Finite_element_method_in_structural_mechanics#Interpolation_or_shape_functions) and its derivatives, the assembly of the discretised stiffness matrix ([@sec:building]), the numerical solution of the system of equations ([@sec:solving]) and the computation of the gradient of the solution (@sec:stress-computation). Yet, the fact that all these a-priori unconnected steps give rise to pretty pictures that resemble reality is still astonishing to me.
+divert(0)
 
-
-Again, take all this information as coming from a fellow that has already taken such a journey from college’s pencil and paper to real engineering cases involving complex numerical calculations. And developing, in the meantime, both an actual working finite-element [back-end](https://www.seamplex.com/fino) and [front-end](https://www.caeplex.com) from scratch.^[The dean of the engineering school I attended used to say “It is not the same to read than to write manuals, and we should aim at writing them.”]
+Again, take all this information as coming from a fellow that has already taken such a journey from college’s pencil and paper to real engineering cases involving complex numerical calculations. And developing, in the meantime, both an actual working finite-element [back-end](https://www.seamplex.com/fino) and [front-end](https://www.caeplex.com) from scratch.dnl^[The dean of the engineering school I attended used to say “It is not the same to read than to write manuals, and we should aim at writing them.”]
 
 
 ### Tips and tricks
@@ -46,8 +52,10 @@ Again, take all this information as coming from a fellow that has already taken
 There are some useful tricks that come handy when trying to solve a mechanical problem. Throughout this text, I will try to tell you some of them.
 
 
-One of the most important ones to use your _imagination_. You will need a lot of imagination to “see” what it is actually going on when analysing an engineering problem. This skill comes from my background in nuclear engineering where I had no choice but to imagine a [positron-electron annihilation](https://en.wikipedia.org/wiki/Electron%E2%80%93positron_annihilation) or an [spontaneous fission](https://en.wikipedia.org/wiki/Spontaneous_fission). But in mechanical engineering, it is likewise important to be able to imagine how the loads “press” one element with the other, how the material reacts depending on its properties, how the nodal displacements generate stresses (both normal and shear), how results converge, etc. And what these results actually mean besides the pretty-coloured figures.^[A former boss once told me “I need the CFD” when I handed in some results. I replied that I did not do computational fluid-dynamics but computed the neutron flux kinetics within a nuclear reactor core. He joked “I know, what I need are the _Colors For Directors_, those pretty-coloured figures along with your actual results to convince the managers.”]
-This journey will definitely need your imagination. We will take a look at equations, numbers, plots, schematics, CAD geometries, 3D\ views, etc. Still, when the theory says “thermal expansion produces normal stresses” you have to picture in your head three little arrows pulling away from the same point in three directions, or whatever mental picture you have about what you understand thermally-induced stresses are. What comes to your mind when someone says that out of the nine elements of the [stress tensor](https://en.wikipedia.org/wiki/Cauchy_stress_tensor) ([@sec:tensor]) there are only six that are independent? Whatever it is, try to practice that kind of graphical thoughts with every new concept. Nevertheless, there will be particular locations of the text where imagination will be most useful. I will bring the subject up then and again.
+One of the most important ones to use your _imagination_. You will need a lot of imagination to “see” what it is actually going on when analysing an engineering problem. This skill comes from my background in nuclear engineering where I had no choice but to imagine a [positron-electron annihilation](https://en.wikipedia.org/wiki/Electron%E2%80%93positron_annihilation) or an [spontaneous fission](https://en.wikipedia.org/wiki/Spontaneous_fission). But in mechanical engineering, it is likewise important to be able to imagine how the loads “press” one element with the other, how the material reacts depending on its properties, how the nodal displacements generate stresses (both normal and shear), how results converge, etc.
+dnl And what these results actually mean besides the pretty-coloured figures.^[A former boss once told me “I need the CFD” when I handed in some results. I replied that I did not do computational fluid-dynamics but computed the neutron flux kinetics within a nuclear reactor core. He joked “I know, what I need are the _Colors For Directors_, those pretty-coloured figures along with your actual results to convince the managers.”]
+
+This journey will definitely need your imagination. We will peek a little bit into equations, numbers, plots, schematics, CAD geometries, 3D\ views, etc. Still, when the theory says “thermal expansion produces normal stresses” you have to picture in your head three little arrows pulling away from the same point in three directions, or whatever mental picture you have about what you understand thermally-induced stresses are. What comes to your mind when someone says that out of the nine elements of the [stress tensor](https://en.wikipedia.org/wiki/Cauchy_stress_tensor) ([@sec:tensor]) there are only six that are independent? Whatever it is, try to practice that kind of graphical thoughts with every new concept. Nevertheless, there will be particular locations of the text where imagination will be most useful. I will bring the subject up then and again.
 
 
 Another heads up is that we will be digging into some mathematics. Probably they would be simple and you would deal with them very easily. But chances are you do not like equations. No problem! Just ignore them for now. Read the text skipping them, it should work as well.
@@ -67,13 +75,13 @@ After further years passed by, engineers (probably the same people that forked s
 
 ![Three-dimensional CAD model for the piping system of\ [@fig:isometric]](real-piping.png){#fig:real-life width=75%}
 
+
 ### Nuclear reactors
 
 In each of the countries that have at least one nuclear power plant there exists a national regulatory body who is responsible for allowing the owner to operate the reactor. These operating licenses are time-limited, with a range that can vary from 25 to 60 years, depending on the design and technology of the reactor. Once expired, the owner might be entitled to an extension, which the regulatory authority can accept provided it can be shown that a certain (and very detailed) set of safety criteria are met. One particular example of requirements is that of fatigue in pipes, especially those that belong to systems that are directly related to the reactor safety.
 
 ### Pressurised pipes
 
-
 How come that pipes are subject to fatigue? Well, on the one hand and without getting into many technical details, the most common nuclear reactor design uses liquid water as coolant and moderator. On the other hand, nuclear power plants cannot by-pass the thermodynamics of the [Carnot cycle](https://en.wikipedia.org/wiki/Carnot_cycle), and in order to maximise the efficiency of the conversion of the energy stored in the uranium nuclei into electricity we need to reach temperatures as high as possible. So, if we want to have liquid water in the core as hot as possible, we need to increase the pressure. The limiting temperature and pressure are given by the [critical point of water](https://en.wikipedia.org/wiki/Critical_point_(thermodynamics)), which is around 374ºC and 22\ MPa. It is therefore expected to have temperatures and pressures near those values in many systems of the plant, especially in the primary circuit and those that directly interact with it, such as pressure and inventory control system, decay power removal system, feedwater supply system, emergency core-cooling system, etc.
 
 Nuclear power plants are not always working at 100% of their maximum power capacity. They need to be maintained and refuelled, they may undergo operational (and some incidental) transients, they might operate at a lower power due to load following conditions, etc. These transient cases involve changes both in temperatures and in pressures that the pipes are subject to, which in turn give rise to changes in the stresses within the pipes. As the transients are postulated to occur cyclically during a number of times during the life-time of the plant (plus its extension period), mechanical fatigue in these piping systems may arise, especially at the interfaces between materials with different thermal expansion coefficients.
@@ -138,7 +146,7 @@ What does this all have to do with mechanical engineering? Well, once we know wh
 
 ### An infinitely-long pressurised pipe {#sec:infinite-pipe}
 
-Let us proceed to our second step, and consider the infinite pipe subject to uniform internal pressure already introduced in\ [@fig:infinite-pipe]. Actually, we are going to solve the mechanical problem on an infinite hollow cylinder, which looks like pipe. This case is usually tackled in college courses, and chances are you already solved it. In fact, the first (and simpler) problem is the “thin cylinder problem.” Then, the “thick cylinder problem” is introduced (the one we solve below), which is slightly more complex. Nevertheless, it has an [analytical solution](https://www.seamplex.com/fino/doc/pipe-linearized/). For the present case, let us consider an infinite pipe (i.e. a hollow cylinder) of internal radius $a$ and external radius $b$ with uniform mechanical properties---Young’s modulus $E$ and Poisson’s ratio $\nu$---subject to an internal uniform pressure $p$.
+Let us proceed to our second step, and consider the infinite pipe subject to uniform internal pressure already introduced in\ [@fig:infinite-pipe]. Actually, we are going to solve the mechanical problem on an infinite hollow cylinder, which looks like pipe. This case is usually tackled in college courses, and chances are you already solved it. In fact, the first (and simpler) problem is the “thin cylinder problem.” Then, the “thick cylinder problem” is introduced (the one we solve below), which is slightly more complex. Nevertheless, it has an analytical solution.^[A detailed analysis of the analytical solution and how results obtained with finite elements converge with respect to the mesh size can be found in the report “On convergence of linearized stresses in an infinite pipe computed using the finite element method.” See [@sec:online].] For the present case, let us consider an infinite pipe (i.e. a hollow cylinder) of internal radius $a$ and external radius $b$ with uniform mechanical properties---Young’s modulus $E$ and Poisson’s ratio $\nu$---subject to an internal uniform pressure $p$.
 
 What follows is more or less what we are taught in school: some equations with a brief explanation of the results. And then we move on to the next subject.
 
@@ -160,7 +168,7 @@ What does this mean? Well, that overall the whole pipe expands a little bit radi
  1. linearly with the pressure, i.e. twice the pressure means twice the displacement, and
  2. inversely proportional to the Young’s Modulus\ $E$ divided by $1+\nu$, i.e. the more resistant the material, the less radial displacements.
 
-That is how an infinite pipe withstands internal pressure.
+That is how an infinite pipe withstands internal pressure. And that is what we are taught in college, which is actually true by the way!
 
 #### Stresses
 
@@ -333,7 +341,7 @@ Remember the main issue of the fatigue analysis in these systems is to analyse w
 
 
 ::::: {#fig:weldolet}
-![Overall view](weldolet-cad1.png){#fig:weldolet-cad width=65%}
+![Overall view](weldolet-cad1-commented.png){#fig:weldolet-cad width=65%}
 
 ![Unstructured tetrahedra-based grid](weldolet-mesh2.png){#fig:weldolet-mesh width=85%}
 
@@ -344,7 +352,7 @@ CAD model of a piping system with a 3/4-inch weldolet-type fork (stainless steel
 ![Location of six SCLs defined to analyse fatigue around a junction.](weldolet-scls-n.png){#fig:weldolet-scls width=75%}
 
 
-On the one hand, a reasonable number of nodes (it is the number of nodes that defines the problem size, not the number of elements as discussed in [@sec:elements-nodes]) in order to get a decent grid is around 200k for each system. On the other hand, solving a couple of dozens of transient heat transfer problems (which we cannot avoid due to the large thermal inertia of the pipes) during a few thousands of seconds over a couple hundred of thousands of nodes might take more time and storage space to hold the results than we might expect.
+On the one hand, a reasonable number of nodes (it is the number of nodes that defines the problem size, not the number of elements) in order to get a decent grid is around 200k for each system. On the other hand, solving a couple of dozens of transient heat transfer problems (which we cannot avoid due to the large thermal inertia of the pipes) during a few thousands of seconds over a couple hundred of thousands of nodes might take more time and storage space to hold the results than we might expect.
 
 There is a wonderful essay by [Isaac Asimov](https://en.wikipedia.org/wiki/Isaac_Asimov) called [“The Relativity of Wrong”](https://en.wikipedia.org/wiki/The_Relativity_of_Wrong) where he introduces the idea that even if something cannot be computed exactly, there are different levels of error. For instance, believing that the Earth is a sphere is less wrong than believing that the Earth is flat, but wrong nonetheless, since it really deviates from a perfect sphere and resembles more an oblate spheroid.
 
@@ -358,7 +366,7 @@ We can then merge this idea by Asimov with an adapted version of the [Saint-Vena
 
 As an example, let us consider the system depicted in\ [@fig:valve] where there is a stainless-carbon steel interface at the discharge of the valves. Instead of solving the transient heat-conduction problem with the internal temperature of the pipes equal to the temperature of the water in the reference transient condition of the power plant and an external condition of natural convection to the ambient temperature in the whole mesh, a reduced model consisting of half of one of the two valves and a small length of the pipes at both the valve inlet and outlet is used. Once the temperature distribution\ $\hat{T}(\mathbf{x},t)$ for each time is obtained in the reduced mesh ([@fig:valve-temp], which has the origin at the centre of the valve), the actual temperature distribution\ $T(\mathbf{x},t)$ is computed by an algebraic generalisation of $\hat{T}(\mathbf{x},t)$ in the full coordinate system. As stated above, those locations which are not covered by the reduced model are generalised with a time-dependent uniform temperature which is the average of the inner and outer temperatures at the inlet and outlet of the reduced mesh.
 
-![Reduced mesh around a valve refined around the interface where the transient heat conduction problem is solved.](valve-temp.png){#fig:valve-temp width=80%}
+![Mesh of half-a-symmetric-valve refined around the interface where the transient heat conduction problem is solved.](valve-temp-commented.png){#fig:valve-temp width=80%}
 
 
 Note that there is no need to have a one-to-one correspondence between the elements from the reduced mesh with the elements from the original one. Actually, the reduced mesh contains first-order elements whilst the former has second-order elements. Also the grid density is different. Nevertheless, the finite-element solver [Fino](https://www.seamplex.com/fino) used to solve both the heat and the mechanical problems, allows to read functions of space and time defined over one mesh and continuously evaluate and use them into another one even if the two grids have different elements, orders or even dimensions. In effect, in the system from\ [@fig:real-life] the material interface is between a orifice plate made in stainless steel that is welded to a carbon-steel pipe ([@fig:real]). The thermal problem can be modelled using a two-dimensional axi-symmetric grid and then generalised to the full three-dimensional mesh using the algebraic manipulation capabilities provided by [Fino](https://www.seamplex.com/fino) (actually by [wasora](https://www.seamplex.com/wasora)) as shown in\ [@fig:real].
@@ -366,16 +374,17 @@ Note that there is no need to have a one-to-one correspondence between the eleme
 dnl ![Full 3D mesh, carbon steel is magenta, stainless steel is green.](real-mesh2.png){#fig:real-mesh2 width=70%}
 dnl ![Temperature distribution for a certain time within the transient computed on a reduced two-dimensional axi-symmetric mesh modelling half the orifice plate and a length of the carbon pipe.](real-temp.png){#fig:real-temp width=34%}
 
-![Temperature distribution for a certain time within the transient computed on a reduced two-dimensional axi-symmetric mesh modelling half the orifice plate and a length of the carbon pipe and Generalisation of the temperature to the full three-dimensional grid.](real-gen.png){#fig:real width=85%}
+![Temperature distribution for a certain time within the transient computed on a reduced two-dimensional axi-symmetric mesh modelling half the orifice plate and a length of the carbon pipe and generalisation of the temperature to the full three-dimensional grid.](real-gen.png){#fig:real width=85%}
 
 
 ### Seismic loads {#sec:seismic}
 
-Before considering the actual mechanical problem that will give us the stress tensor at the SCLs, and besides needing to solve the transient thermal problem to get the temperature distributions, we need to address the loads that arise due to a postulated earthquake during a certain part of the operational transients. The full computation of a mechanical transient problem using the earthquake time-dependent displacements is off the table for two reasons. First, because again the computation would take more time than we might have to deliver the report. And secondly and more importantly, because civil engineers do not compute earthquakes in the time domain but in the frequency domain using the [response spectrum method](https://en.wikipedia.org/wiki/Response_spectrum), provided the case is linear (see\ [@sec:linearity]). Time to revisit our [Fourier transform](https://en.wikipedia.org/wiki/Fourier_transform) exercises from undergraduate math courses.
 
+Every nuclear power plant is designed to withstand earthquakes. Of course, not all plants need the same level of reinforcements. Those built in large quiet plains will be, seismically speaking, cheaper than those located in geologically active zones. Keep in mind that all the 54 Japanese nuclear power plants did structurally resist the 2011 earthquake, and all of the reactors were safely shut down. What actually happened in [Fukushima](http://www.world-nuclear.org/information-library/safety-and-security/safety-of-plants/fukushima-accident.aspx) is that one hour after the main shake, a 14-metre tsunami splashed on the coast, jumping over the 9-metre defences and flooding the emergency [Diesel generators](https://en.wikipedia.org/wiki/Diesel_generator) that provided power to the pumps in charge of removing the remaining [decay power](https://en.wikipedia.org/wiki/Decay_heat) from the already-stopped [reactor core](https://en.wikipedia.org/wiki/Nuclear_reactor_core).
+
+divert(-1)
 #### Earthquake spectra
 
-Every nuclear power plant is designed to withstand earthquakes. Of course, not all plants need the same level of reinforcements. Those built in large quiet plains will be, seismically speaking, cheaper than those located in geologically active zones. Keep in mind that all the 54 Japanese nuclear power plants did structurally resist the 2011 earthquake, and all of the reactors were safely shut down. What actually happened in [Fukushima](http://www.world-nuclear.org/information-library/safety-and-security/safety-of-plants/fukushima-accident.aspx) is that one hour after the main shake, a 14-metre tsunami splashed on the coast, jumping over the 9-metre defences and flooding the emergency [Diesel generators](https://en.wikipedia.org/wiki/Diesel_generator) that provided power to the pumps in charge of removing the remaining [decay power](https://en.wikipedia.org/wiki/Decay_heat) from the already-stopped [reactor core](https://en.wikipedia.org/wiki/Nuclear_reactor_core).
 
 Back to our case study, the point is that each site where nuclear power plants are built must have a geological study where a postulated design-basis earthquake is to be defined. In other words, a theoretical earthquake which the plant ought to withstand needs to be specified. How? By giving a set of three spectra (one for each coordinate direction) giving acceleration as a function of the frequency for each level of the building. That is to say, once the earthquake hits the power plant, depending on soil-structure interactions the energy will shake the building foundations in a way that depends on the characteristics of the earthquake, the soil and the concrete structure. Afterwards, the way the oscillations travel upward and shake each of the mechanical components erected in each floor level depends on the design of the civil structure in a way which is fully determined by the floor response spectra like the ones depicted in\ [@fig:spectrum].
 
@@ -399,8 +408,13 @@ A real continuous solid has infinite modes of oscillation. A discretised one (us
 
 ![Effective accumulated mass fraction up to the $i$-th mode of oscillation.](acumulada.emf){#fig:acumulada}
 
+
+These first modes, shown in\ [@fig:modes], that take up most of the energy are then used to take into account the earthquake load. There are several ways of performing this computation, but the ASME\ III code states that the method known as SRSS (for Square Root of Sum of Squares) can be used. This method mixes the eigenvectors with the floor response spectra through the eigenvalues and gives a spatial (actually nodal) distribution of three accelerations (one for each direction) that, when multiplied by the density of the material, give a vector of a distributed force (in units of Newton per cubic millimetre for example) which is statically equivalent to the load coming from the postulated earthquake.
+divert(0)
+
+
 ::::: {#fig:modes}
-![$f_1 \approx 30$\ [Hz](https://www.seamplex.com/docs/nafems4/mode1.webm)](mode1.png){width=48%}\  
+![$f_1 \approx 30$\ [Hz](https://www.seamplex.com/docs/nafems4/mode1.webm)](mode1-commented.png){width=48%}\  
 ![$f_2 \approx 60$\ [Hz](https://www.seamplex.com/docs/nafems4/mode2.webm)](mode2.png){width=48%}
 
 ![$f_3 \approx 70$\ [Hz](https://www.seamplex.com/docs/nafems4/mode3.webm)](mode3.png){width=48%}\  
@@ -409,12 +423,9 @@ A real continuous solid has infinite modes of oscillation. A discretised one (us
 ![$f_5 \approx 100$\ [Hz](https://www.seamplex.com/docs/nafems4/mode5.webm)](mode5.png){width=48%}\  
 ![$f_6 \approx 130$\ [Hz](https://www.seamplex.com/docs/nafems4/mode6.webm)](mode6.png){width=48%}
 
-First six natural oscillation modes for a piping section between axial supports. Blue-red color scale shows magnitude of deformation, as compared to the original geometry (grey). Animations can be seen at links like <https://www.seamplex.com/docs/nafems4/mode1.webm>.
+First six natural oscillation modes for a piping section between axial supports. Blue-red colour scale shows magnitude of deformation, as compared to the original geometry (grey). Full animations available online ([@sec:online]).
 :::::
 
-
-These first modes, shown in\ [@fig:modes], that take up most of the energy are then used to take into account the earthquake load. There are several ways of performing this computation, but the ASME\ III code states that the method known as SRSS (for Square Root of Sum of Squares) can be used. This method mixes the eigenvectors with the floor response spectra through the eigenvalues and gives a spatial (actually nodal) distribution of three accelerations (one for each direction) that, when multiplied by the density of the material, give a vector of a distributed force (in units of Newton per cubic millimetre for example) which is statically equivalent to the load coming from the postulated earthquake.
-
 ::::: {#fig:acceleration}
 ![$a_x$](ax.png){width=80%}
 
@@ -422,10 +433,11 @@ These first modes, shown in\ [@fig:modes], that take up most of the energy are t
 
 ![$a_z$](az.png){width=80%}
 
-The equivalent accelerations for the piping section of [@fig:modes] for the spectra of\ [@fig:spectrum]
+Equivalent accelerations for a certain piping section
 :::::
 
-The ASME code says that these accelerations (depicted in [@fig:acceleration]) ought to be applied twice: once with the original sign and once with all the elements with the opposite sign, as SRSS “looses signs.” The application of each of these equivalent loads should last two seconds in the original time domain.
+Since the computation of the loads that a certain earthquake gives rise to would add a significant amount of complexity to this already-complex case study, the details are skipped. In any case, we need to compute the first natural modes of oscillation of the piping section under study ([@fig:modes]) and then, after some cumbersome maths, obtain an statically-equivalent load distribution. In effect, [@fig:acceleration] shows a sample distribution of acceleration distribution within a certain piping system which, when multiplied by the density, give a load distribution which is statically equivalent to the dynamic solicitations created by the earthquake.
+The ASME code says that these accelerations ought to be applied twice: once with the original sign and once with all the elements with the opposite sign, as SRSS “looses signs.” The application of each of these equivalent loads should last two seconds in the original time domain.
 
 
 ### Linearity (not yet linearisation) {#sec:linearity}
@@ -436,17 +448,17 @@ Let us jump out of our nuclear piping problem and step back into the general fin
 
 \medskip
 
-Time for both of us to make an experiment. Grab your favourite FEM program for the first time (remember mine is [CAEplex](https://caeplex.com)) and create a 1mm $\times$ 1mm $\times$ 1mm cube. Set any values for the Young’s Modulus and Poisson ratio as you want. I chose\ $E=200$\ GPa and\ $\nu=0.28$. Restrict the three faces pointing to the negative axes to their planes, i.e.
+Time for both of us to make an experiment. Grab your favourite FEM program for the first time (remember mine is [CAEplex](https://caeplex.com), which can also be accessed through [Onshape](https://www.onshape.com/cad-blog/partner-spotlight-using-caeplex-for-finite-element-analysis-in-onshape)) and create a 1mm $\times$ 1mm $\times$ 1mm cube. Set any values for the Young’s Modulus and Poisson ratio as you want. I chose\ $E=200$\ GPa and\ $\nu=0.28$. Restrict the three faces pointing to the negative axes to their planes, i.e.
 
  * in face “left” ($x<0$), set null displacement in the $x$ direction ($u=0$),
  * in face “front” ($y<0$), set null displacement in the $y$ direction ($v=0$),
  * in face “bottom” ($z<0$), set null displacement in the $z$ direction ($w=0$).
 
-Now we are going to create and compare three load cases:
+Now we are going to create and compare three load cases:^[The provided links lead to FEA cases which are fully usable directly from the web browser, i.e. “finite elements in the cloud.”]
 
- a. Pure normal loads (<https://caeplex.com/p/d8f>)
+ a. Pure normal loads (<https://caeplex.com/p/d8fe>)
  b. Pure shear loads (<https://caeplex.com/p/b494>)
- c. The combination of A & B (<https://caeplex.com/p/989>)
+ c. The combination of A & B (<https://caeplex.com/p/9899>)
 
 The loads in each cases are applied to the three remaining faces, namely “right” ($x>0$), “back” ($y>0$) and “top,” ($z>0$). Their magnitude in Newtons are:
 
@@ -477,8 +489,9 @@ $$
 
 
 ::::: {#fig:cube}
-![Case B, pure-shear loads](cube-shear.png){#fig:cube-shear width=48%}\  
-![Case C, normal plus shear loads](cube-full.png){#fig:cube-full width=48%}
+![Case B, pure-shear loads](cube-shear.png){#fig:cube-shear width=65%}
+
+![Case C, normal plus shear loads](cube-full.png){#fig:cube-full width=65%}
 
 Spatial distribution of principal stress\ 3 for cases\ B and\ C. If linearity applied, case\ C would be equal to case\ B plus a constant
 :::::
@@ -623,14 +636,14 @@ The first thing that has to be said is that, as with any interesting problem, th
  10. etc.
 
 ::::: {#fig:quarter}
-![Structured second-order incomplete hexahedra](quarter-struct.png){#fig:cube-struct width=48%}\  
-![Unstructured second-order tetrahedra](quarter-caeplex.png){#fig:quarter-caeplex width=48%}
+![Structured second-order hexahedra](quarter-struct.png){#fig:cube-struct width=50%}\  
+![Unstructured second-order tetrahedra](quarter-caeplex.png){#fig:quarter-caeplex width=50%}
 
 Two of the hundreds of different ways the infinite pressurised pipe can be solved using FEM. The axial displacement at the ends is set to zero, leading to a [plane-strain](https://en.wikipedia.org/wiki/Plane_stress#Plane_strain_(strain_matrix)) condition
 :::::
 
 
-You can get both the exponential nature of each added bullet and how easily we can add new further choices to solve a FEM problem. And each of these choices will reveal you something about the nature of either the mechanical problem or the numerical solution. It is not possible to teach any possible lesson from every outcome in college, so you will have to learn them by yourself getting your hands at them. I have already tried to address the particular case of the infinite pipe in a [recent report](https://www.seamplex.com/fino/doc/pipe-linearized/)^[<https://www.seamplex.com/fino/doc/pipe-linearized/>] that is worth reading before carrying on with this article.
+You can get both the exponential nature of each added bullet and how easily we can add new further choices to solve a FEM problem. And each of these choices will reveal you something about the nature of either the mechanical problem or the numerical solution. It is not possible to teach any possible lesson from every outcome in college, so you will have to learn them by yourself getting your hands at them. I have already tried to address the particular case of the infinite pipe in a [recent report](https://www.seamplex.com/fino/doc/pipe-linearized/)^[The aforementioned “On convergence of linearized stresses in an infinite pipe computed using the finite element method.”] that is worth reading before carrying on with this article.
 The main conclusions of the report are:
 
 dnl * Engineering problems ought not to be solved using black-boxes (i.e. privative software whose source code is not freely available)---more on the subject below in\ [@sec:two-materials].
@@ -640,12 +653,12 @@ dnl * There are no shear stresses, so these three stresses are also the principa
 dnl * Analytical expressions for the membrane and membrane plus bending stresses along any radial SCL can be obtained.
 dnl * The spatial domain can be discretised using linear or higher-order elements. In particular first and second-order elements have been used in the report.
  * For the same number of elements, second-order grids need more nodes than linear ones, although they can better represent curved geometries.
- * The discretised problem size depends on the number of nodes and not on the number of elements---more on the subject below in\ [@sec:elements-nodes].
+dnl * The discretised problem size depends on the number of nodes and not on the number of elements---more on the subject below in\ [@sec:elements-nodes].
 dnl * The finite-element results for the displacements are similar to the analytical solution, with second-order grids giving better results for the same number of elements (this is expected as they involved far more nodes).
  * The three stress distributions computed with the finite-element give far more reasonable results for the second-order case than for the first-order grid. This is qualitatively explained by the fact that  first-order tetrahedra have uniform derivatives and such the elements located in both the external and external faces represent the stresses not at the actual faces but at the barycentre of the elements.
  * Membrane stresses converge well for both the first and second-order cases because they represent a zeroth-order moment of the stress distribution and the excess and defect errors committed by the internal and external elements approximately cancel out.
  * Membrane plus bending stresses converge very poorly with linear elements because the excess and defect errors do not cancel out because it is a first-order moment of the stress distribution.
- * The computational effort to solve a given problem, namely the CPU time and the needed RAM ([@fig:error-vs-cpu]) depend non-linearly on various factors, but the most important one is the problem size which is three times the number of nodes in the grid---more on the subject below in\ [@sec:elements-nodes].
+ * The computational effort to solve a given problem, namely the CPU time and the needed RAM ([@fig:error-vs-cpu]) depend non-linearly on various factors, but the most important one is the problem size which is three times the number of nodes in the grid dnl---more on the subject below in\ [@sec:elements-nodes].
  * The error with respect to the analytical solutions as a function of the CPU time needed to compute the membrane stress is similar for both first and second-order grids. But for the computation of the membrane plus bending stress ([@fig:error-MB-vs-cpu]), first-order grids give very poor results compared to second-order grids for the same CPU time.
 
 ::::: {#fig:error-vs-cpu}
@@ -657,6 +670,7 @@ Error in the computation of the linearised stresses vs. CPU time needed to solve
 
 An additional note should be added. The FEM solution, which not only gives the nodal displacements but also a method to interpolate these values inside the elements, does not fully satisfy the original equilibrium equations at every point (i.e. the strong formulation). It is an approximation to the solution of the [weak formulation](https://en.wikipedia.org/wiki/Weak_formulation) that is close (measured in the vector space spanned by the [shape functions](https://www.quora.com/What-is-a-shape-function-in-FEM)) to the real solution. Mechanically, this means that the FEM solution leads only to nodal equilibrium but not point-wise equilibrium.
 
+divert(-1)
 ### Elements, nodes and CPU {#sec:elements-nodes}
 
 dnl The last two bullets above lead to an issue that has come many times when discussing the issue of convergence with respect to the mesh size with other colleagues. There apparently exists a common misunderstanding that the number of elements is the main parameter that defines how complex a FEM model is. This is strange, because even in college we are taught that the most important parameter is the _size_ of the stiffness matrix, which is three times (for 3D problems with the [displacement-based formulation](http://web.mit.edu/kjb/www/Books/FEP_2nd_Edition_4th_Printing.pdf) formulation) the number of _nodes_.
@@ -732,7 +746,7 @@ Now proceed to picturing the general three-dimensional cases with unstructured t
 Detailed mathematics show that the location where the derivatives of the interpolated displacements are closer to the real (i.e. the analytical ones in problems that have it) solution are the elements’ [Gauss points](https://en.wikipedia.org/wiki/Gaussian_quadrature). Even better, the material properties at these points are continuous (they are usually uniform but they can depend on temperature for example) because, unless we are using weird elements, there are no material interfaces inside elements. But how to manage a set of stresses given at the Gauss points instead of at the nodes? Should we use one mesh for the input and another one for the output? What happens when we need to know the stresses on a surface and not just in the bulk of the solid? There are still no one-size-fits-all answers. There is a very interesting [blog post](https://nickjstevens.netlify.com/post/2019/stress-singularities/) by Nick Stevens that addresses the issue of stresses computed at sharp corners. What does your favourite FEM program do with such a case?
 
 In any case, this step takes a non-negligible amount of time. The most-common approach, i.e. the node-averaging method is driven mainly by the number of nodes of course. So all-in-all, these are the reasons to use the number of nodes instead of the numbers of elements as a basic parameter to measure the complexity of a FEM problem.
-
+divert(0)
 
 ## Adding complexity: the truth is out there
 
@@ -861,7 +875,7 @@ Do you now see the added value of training throw-ins with watermelons? We might
  * using distributed forces from earthquakes,
  * etc.
 
-Most of the time at college we would barely do what is needed to be approved in one course and move on to the next one. If you have the time and consider a career related to finite-element analysis, please do not. Clone the repository^[<https://github.com/seamplex/tee>] with the input files for [Fino](https://www.seamplex.com/fino) and start playing. If you are stuck, do not hesitate asking for help in [wasora’s mailing list](https://www.seamplex.com/lists.html).
+Most of the time at college we would barely do what is needed to be approved in one course and move on to the next one. If you have the time and consider a career related to finite-element analysis, please do not. Clone the repository ([@sec:online]) with the input files for [Fino](https://www.seamplex.com/fino) and start playing. If you are stuck, do not hesitate asking for help in [wasora’s mailing list](https://www.seamplex.com/lists.html).
 
 One further detail: it is always a sane check to try to explain the numerical results based on physical reasoning (i.e. “with your fingers”) as we did two paragraphs above. Most of the time you will be solving problems whilst already knowing what the result would (or ought to) be.
 
@@ -872,7 +886,7 @@ A fellow mechanical engineer who went to the same high school I did, who went to
 
 So here we are again with the case study where we have to compute the linearised stresses at certain SCLs located near the interface between two different kinds of steels during operational and incidental transients of the plant. The first part is then to “bake” the pipes, modelling a thermal transient with time-dependent temperature or convection boundary conditions (it depends on the system). This steps gives a time and space-dependent temperature\ $T(x,y,z,t)$.
 
-Then we proceed to “shake” the pipes, i.e. to compute the natural frequencies and associated oscillations modes. Employing the floor response spectra and combining the individual contributions with the SRSS method discussed in section\ [@sec:seismic], we obtain a distributed load vector\ $\mathbf{f}(x,y,z)$ which is statically equivalent to the design earthquake.
+We then proceed to “shake” the pipes. That is to say, we obtain a distributed load vector\ $\mathbf{f}(x,y,z)$ which is statically equivalent to the design earthquake.
 
 Finally we attempt to “break” the pipes successively solving many steady-state elastic problems for different times\ $t$ of the operational transient. Some remarks about this step:
 
@@ -908,7 +922,11 @@ If you notice, the complex dependence of the stiffness matrix\ $K$ can be reduce
 
 $$ K(\mathbf{x}) \cdot \mathbf{u}(\mathbf{x}) = \mathbf{b}(\mathbf{x})$$
 
-And this last equation is linear in\ $\mathbf{u}$. In effect, the discretisation step means to integrate over\ $\mathbf{x}$. As\ $K$, $\mathbf{u}$ and\ $\mathbf{b}$ depend only on\ $\mathbf{x}$, then after integration one gets just numbers with the matrix representation of\ [@eq:kub]. Again, you can either trust me, ask a teacher or go through with the maths (in increasing order of recommendation).
+And this last expression is linear in\ $\mathbf{u}$! In effect, the discretisation step means to integrate over\ $\mathbf{x}$. As\ $K$, $\mathbf{u}$ and\ $\mathbf{b}$ depend only on\ $\mathbf{x}$, then after integration one gets just numbers inside $K$ and $\mathbf{b}$. Again, you can either, in increasing order of recommendation:
+
+ 1. trust me,
+ 2. ask a teacher, or
+ 3. go through with the maths.
 
 
 ::::: {#fig:case}
@@ -943,7 +961,7 @@ Temperature distribution for a certain instant of the transient, computed in the
 ![$f_1 \approx 930$\ [Hz](https://www.seamplex.com/docs/nafems4/case-mode9.webm)](case-mode9.png){#fig:case-mode-9 width=30%}
 
 First nine natural modes of oscillation of the piping system subject to the boundary conditions the supports provide.
-Animations can be seen at links like <https://www.seamplex.com/docs/nafems4/case-mode1.webm>.
+Animations are available online ([@sec:online]).
 :::::
 
 dnl ![Floor response spectrum.](case-spectrum.emf){#fig:case-spectrum width=70%}
@@ -955,15 +973,15 @@ To recapitulate, the figures in this section show some partial non-dimensional r
  3. defining the number and locations of the SCLs ([@fig:case-scls] or [@fig:weldolet-scls;@fig:valve])
  4. computing a heat conduction (bake) transient problem with temperatures as a function of time from the operational transient in a simple domain using temperature-dependent thermal conduction coefficients from ASME\ II\ part D ([@fig:case-temp] or [@fig:valve-temp])
  5. generalising the temperature distribution as a function of time to the general domain ([@fig:case-temp2] or [@fig:real])
- 6. performing a modal analysis (shake) on the main domain to obtain the main oscillation frequencies and modes ([@fig:case-mode] or [@fig:modes])
- 7. using the floor response spectra ([@fig:spectrum]) and the SRSS method to obtain a distributed force statically-equivalent to the earthquake load ([@fig:case-acceleration] or [@fig:acceleration])
+ 6. performing a modal analysis (shake) on the main domain to obtain the main oscillation frequencies and modes ([@fig:case-mode] dnl or [@fig:modes])
+ 7. obtaining a distributed force statically-equivalent to the earthquake load (shake, [@fig:case-acceleration] or [@fig:acceleration])
  8. successively solving the linear elastic problem for different times using the generalised temperature distribution taking into account
      a. the dependence of the Young’s Modulus\ $E$ and the thermal expansion coefficient\ $\alpha$ with temperature,
      b. the thermal expansion effect itself
      c. the instantaneous pressure exerted in the internal faces of the pipes at the time\ $t$ according to the definition of the operational transient
      d. the restriction of the degrees of freedom of those faces, lines or points that correspond to mechanical supports located both within and at the ends of the CAD model
      e. the earthquake load, which according to ASME should be present only during four seconds of the transient: two seconds with one sign and the other two seconds with the opposite sign. This period should be selected to coincide with the instant of highest mechanical stress (conservative computation)
- 9. computing the linearised stresses (membrane and membrane plus bending) at the SCLs combining them as Principal\ 1,  Principal\ 2, Principal\ 3 and Tresca
+ 9. computing the linearised stresses (membrane and membrane plus bending) at the SCLs combining them as Principal\ 1,  Principal\ 2, Principal\ 3 and Tresca^[A deep technical note should be added for discussing the validity of linearising a stress invariant as the principall stresses individually. It should be enough for the present study to keep in mind that in [@sec:in-air] it is the difference of these linearised principal stresses (i.e. Tresca-like) are used to compute the stress amplitude of the periodic load for fatigue assessment.]
  10. juxtaposing these linearised stresses for each time of the transient and for each transient so as to obtain a single time-history of stresses including all the operational and/or incidental transients under study, which is what stress-based fatigue assessment needs (recall\ [@sec:fatigue] and go on to\ [@sec:usage]).
 
 A pretty nice list of steps, which definitely I would not have been able to tackle when I was in college. Would you?
@@ -1009,7 +1027,8 @@ When\ $\text{CUF} < 1$, then the part under analysis can withstand the proposed
 This cryptic paragraph can be better explained by using a clearer example. To avoid using actual sensitive data from a real power plant, let us use the same test case used by both the [US Nuclear Regulatory Commission](https://en.wikipedia.org/wiki/Nuclear_Regulatory_Commission) (in its report NUREG/CR-6909) and the [Electric Power Institute](https://en.wikipedia.org/wiki/Electric_Power_Research_Institute) (report 1025823) called “EAF (Environmentally-Assisted Fatigue) Sample Problem 2-Rev.\ 2 (10/21/2011)”.
 
 
-![A low-alloy steel vessel nozzle (blue) welded to a stainless steel pipe (grey)](axi-inches-3d.png){#fig:axi-inches-3d width=25%}
+![A low-alloy steel vessel nozzle (blue) welded to a stainless steel pipe (grey)](axi-inches-3d.png){#fig:axi-inches-3d width=35%}
+
 
 
 It consists of a typical vessel nozzle with attached piping as illustrated in\ [@fig:axi-inches-3d]. These components are subject to four transients\ $k=1,2,3,4$ that give rise to linearised stress histories (slightly modified according to NB-3216.2) which are given as individual stress values juxtaposed so as to span a time range of about 36,000 seconds ([@fig:nureg1]). As the time steps is not constant, each stress value has an integer index\ $i$ that uniquely identifies it:
@@ -1021,8 +1040,24 @@ It consists of a typical vessel nozzle with attached piping as illustrated in\ [
 |   3   |  6450--9970    |  960--1595  |         20       |
 |   4   |  9970--35971   |  1596--2215 |        100       |
 
-A design-basis earthquake was assumed to occur briefly during one second (sic) at around\ $t=3470$\ seconds, and it is assigned a number of cycles\ $n_e=5$. The detailed stress history for one of the SCLs including both the principal and lineariased stresses, which are already offset following NB-3216.2 so as to have a maximum stress equal to zero, can be found as an appendix in NRC's NUREG/CR-6909 report, or in the repository with the scripts I prepared for you to play with this problem.^[<https://github.com/seamplex/cufen>]
 
+A design-basis earthquake was assumed to occur briefly during one second (sic) at around\ $t=3470$\ seconds, and it is assigned a number of cycles\ $n_e=5$. The detailed stress history for one of the SCLs including both the principal and lineariased stresses, which are already offset following NB-3216.2 so as to have a maximum stress equal to zero, can be found as an appendix in NRC's NUREG/CR-6909 report, or in the repository with the scripts I prepared for you to play with this problem.
+
+To compute the global usage factor, we first need to find all the combinations of local extrema pairs and then sort them in decreasing order of stress difference. For example, the largest stress amplitude is found between $i=447$ and $i=694$. The second one is 447--699. Then 699--1020, 699--899, etc. For each of these pairs, defined by the indexes\ $i_{1,j}$ and $i_{2,j}$, a partial usage factor\ $U_j$ should computed. The stress amplitude\ $S_{\text{alt},j}$ which should be used to enter into the $S$-$N$ curve is
+
+\pagebreak
+
+
+$$
+S_{\text{alt},j} = \frac{1}{2} \cdot k_{e,j} \cdot \left| MB^\prime_{i_{1,j}} - MB^\prime_{i_{2,j}} \right| \cdot \frac{E_\text{SN}}{E(T_{\text{max}_j})}
+$$
+
+\noindent where $k_e$ is a plastic correction factor for large loads (NB-3228.5), $E_\text{SN}$ is the Young’s Modulus used to create the $S$-$N$ curve (provided in the ASME fatigue curves) and\ $E(T_{\text{max}_j})$ is the material’s Young’s Modulus at the maximum temperature within the\ $j$-th interval.
+
+
+We now need to comply with ASME’s obscure note about the number of cycles to assign a proper value of\ $n_j$. Back to the largest pair 447-694, we see that 447 belongs to transient\ #1 which has assigned 20 cycles and 694 belongs to the earthquake with 5 cycles. Therefore $n_1=5$, and the cycles associated to each index are decreased in five. So $i=694$ disappears and the number of cycles associated to $i=447$ are decreased from 20 to 15. The second largest pair is now 447-699, with 15 (because we just spent 5 in the first pair) and 50 cycles respectively. Now $n_2=15$, point 447 disappears and 699 remains with 35 cycles. The next pair is 699-1020, with number of cycles 35 and 20 so $n_3=20$, point 1020 disappears and 699 remains with 15 cycles. And so on, down to the last pair.
+
+Why all these details? Not because I want to teach you how to perform fatigue evaluations just reading this section without resorting to ASME, fatigue books and even other colleagues. It is to show that even though these computation can be made “by hand” (i.e. using a calculator or, God forbids, a spreadsheet) when having to evaluate a few SCLs within several piping systems it is far (and I mean really far) better to automate all these steps by writing a set of scripts. Not only will the time needed to process the information be reduced, but also the introduction of human errors will be minimised and repeatability of results will be assured---especially if working under a [distributed version control](https://en.wikipedia.org/wiki/Distributed_version_control) system such as [Git](https://en.wikipedia.org/wiki/Git). This is true in general, so here is another tip: learn to write scripts to post-process your FEM results (you will need to use a script-friendly FEM program) and you will gain considerable margins regarding time and quality.
 
 ::::: {#fig:nureg}
 ![Full time range of the juxtaposed transients spanning $\approx$ 36,000 seconds](nureg1.emf){#fig:nureg1 width=70%}
@@ -1034,14 +1069,6 @@ A design-basis earthquake was assumed to occur briefly during one second (sic) a
 Time history of the linearised stress\ $\text{MB}_{31}$ corresponding to the example problem from NRC and EPRI. The indexes\ $i$ of the extrema are shown in green (minimums) and red (maximums)
 :::::
 
-
-To compute the global usage factor, we first need to find all the combinations of local extrema pairs and then sort them in decreasing order of stress difference. For example, the largest stress amplitude is found between $i=447$ and $i=694$. The second one is 447--699. Then 699--1020, 699--899, etc. For each of these pairs, defined by the indexes\ $i_{1,j}$ and $i_{2,j}$, a partial usage factor\ $U_j$ should computed. The stress amplitude\ $S_{\text{alt},j}$ which should be used to enter into the $S$-$N$ curve is
-
-$$
-S_{\text{alt},j} = \frac{1}{2} \cdot k_{e,j} \cdot \left| MB^\prime_{i_{1,j}} - MB^\prime_{i_{2,j}} \right| \cdot \frac{E_\text{SN}}{E(T_{\text{max}_j})}
-$$
-where $k_e$ is a plastic correction factor for large loads (NB-3228.5), $E_\text{SN}$ is the Young’s Modulus used to create the $S$-$N$ curve (provided in the ASME fatigue curves) and\ $E(T_{\text{max}_j})$ is the material’s Young’s Modulus at the maximum temperature within the\ $j$-th interval.
-
 ::::: {#fig:cuf}
 ![Reference from NUREG/CR6909](cuf-nrc.png){#fig:cuf-nrc width=100%}
 
@@ -1050,9 +1077,6 @@ where $k_e$ is a plastic correction factor for large loads (NB-3228.5), $E_\text
 Tables of individual usage factors for the NRC/EPRI “EAF Sample Problem 2-Rev.\ 2 (10/21/2011).” One table is taken from a document issued by almost-a-billion-dollar-year-budget government agency from the most powerful country in the world and the other one is from a third-world engineering startup. Guess which is which.
 :::::
 
-We now need to comply with ASME’s obscure note about the number of cycles to assign a proper value of\ $n_j$. Back to the largest pair 447-694, we see that 447 belongs to transient\ #1 which has assigned 20 cycles and 694 belongs to the earthquake with 5 cycles. Therefore\ $n_1=5$, and the cycles associated to each index are decreased in five. So $i=694$ disappears and the number of cycles associated to $i=447$ are decreased from 20 to 15. The second largest pair is now 447-699, with 15 (because we just spent 5 in the first pair) and 50 cycles respectively. Now $n_2=15$, point\ 447 disappears and 699 remains with 35 cycles. The next pair is 699-1020, with number of cycles 35 and 20 so $n_3=20$, point 1020 disappears and 699 remains with 15 cycles. And so on, down to the last pair.
-
-Why all these details? Not because I want to teach you how to perform fatigue evaluations just reading this section without resorting to ASME, fatigue books and even other colleagues. It is to show that even though these computation can be made “by hand” (i.e. using a calculator or, God forbids, a spreadsheet) when having to evaluate a few SCLs within several piping systems it is far (and I mean really far) better to automate all these steps by writing a set of scripts. Not only will the time needed to process the information be reduced, but also the introduction of human errors will be minimised and repeatability of results will be assured---especially if working under a [distributed version control](https://en.wikipedia.org/wiki/Distributed_version_control) system such as [Git](https://en.wikipedia.org/wiki/Git). This is true in general, so here is another tip: learn to write scripts to post-process your FEM results (you will need to use a script-friendly FEM program) and you will gain considerable margins regarding time and quality.
 
 ### In water (NRC’s extension) {#sec:in-water}
 
@@ -1096,7 +1120,7 @@ Tables of individual environmental correction and usage factors for the NRC/EPRI
 
 ## Conclusions
 
-Back in college, we all learned how to solve engineering problems. And once we graduated, we felt we could solve and fix the world (if you did not graduate yet, you will feel it shortly). But there is a real gap between the equations written in chalk on a blackboard (now probably in the form of beamer slide presentations) and actual real-life engineering problems. This chapter introduces a real case from the nuclear industry and starts by idealising the structure such that it has a known analytical solution that can be found in textbooks.  Additional realism is added in stages allowing the engineer to develop an understanding of the more complex physics and a faith in the veracity of the finite-element results where theoretical solutions are not available. Even more, a brief insight into the world of evaluation of stress-life fatigue using such results further illustrates the complexities of real-life engineering analysis. Here is a list of the tips that arose throughout the text:
+Back in college, we all learned how to solve engineering problems. And once we graduated, we felt we could solve and fix the world (if you did not graduate yet, you will feel it shortly). But there is a real gap between the equations written in chalk on a blackboard (now probably in the form of beamer slide presentations) and actual real-life engineering problems. This chapter introduces a real case from the nuclear industry and starts by idealising the structure such that it has a known analytical solution that can be found in textbooks.  Additional realism is added in stages allowing the engineer to develop an understanding of the more complex physics and a faith in the veracity of the finite-element results where theoretical solutions are not available. Even more, a brief insight into the world of evaluation of stress-life fatigue using such results further illustrates the complexities of real-life engineering analysis, even though the presented case was simplified for the sake of clearness. Here is a list of the tips that arose throughout the text:
 
  * use and exercise your imagination
  * practise math
@@ -1105,29 +1129,28 @@ Back in college, we all learned how to solve engineering problems. And once we g
  * remember there are other methods beside finite elements
  * take into account that even within the finite element method, there is a wide variety of complexity in the problems that can be solved
 dnl * follow the “five whys rule” before compute anything, probably you do not need to
- * use engineering judgment and make sure understand Asimov’s [“wronger than wrong”](https://en.wikipedia.org/wiki/Wronger_than_wrong) concept
+dnl * use engineering judgement and make sure understand Asimov’s [“wronger than wrong”](https://en.wikipedia.org/wiki/Wronger_than_wrong) concept
  * play with your favourite FEM solver (mine is [CAEplex](https://caeplex.com)) solving simple cases, trying to predict the results and picturing the stress tensor and its eigenvectors in your imagination
  * prove (with pencil and paper) that even though the principal stresses are not linear with respect to summation, they are linear with respect to multiplication 
  * grab any stress distribution from any of your FEM projects, compute the iso-stress curves and the draw normal lines to them to get acquainted with SCLs
  * search online for “stress linearisation” (or “linearization” if you want) and then get a copy of ASME\ VIII Div\ 2 Annex 5-A
  * take into account that FEM solutions lead only to nodal equilibrium but not point-wise equilibrium
- * measure the time needed to generate grids of different sizes and kinds with your favourite mesher
- * learn this by heart: the complexity of a FEM problem is given mainly by the number of _nodes_, not by the number of elements
+dnl * measure the time needed to generate grids of different sizes and kinds with your favourite mesher
+dnl * learn this by heart: the complexity of a FEM problem is given mainly by the number of _nodes_, not by the number of elements
  * remember that welded materials with different thermal expansion coefficients may lead to fatigue under cyclic temperature changes
 dnl * if you have time, try to get out of your comfort zone and do more than what others expect from you (like parametric computations)
- * clone the [parametric tee repository](https://github.com/seamplex/tee), understand how the figures from\ [@sec:parametric] were built and expand them to cover “we might go on...” bullets
  * try to find an explanation of the results obtained, just like we did when we explained the parametric curves from\ [@fig:tee-MB ] with two opposing effects which were equal in magnitude around $d_b=5$\ inches
- * think thermal-mechanical plus earthquakes as “bake, break and shake” problems
+ * think thermal-mechanical plus earthquakes cases as “bake, break and shake” problems
  * understand why the elastic problem of the case study is still linear after all
  * keep in mind that finite-elements are a means to get an engineering solution, not an end by themselves
- * learn to write scripts to post-process FEM results (from a script-friendly FEM program)
+ * learn to write scripts to post-process FEM results (from an script-friendly open-source FEM program)
  * work under a [distributed version control](https://en.wikipedia.org/wiki/Distributed_version_control) system such as [Git](https://en.wikipedia.org/wiki/Git), even when just editing input files or writing reports
-dnl * clone the [environmental fatigue sample problem repository](https://bitbucket.org/seamplex/cufen) and obtain a nicely-formatted table with the results of the “EAF Sample Problem 2-Rev.\ 2 (10/21/2011)” from\ [@sec:in-air;@sec:in-water].
  * do not write ambiguous reports by replacing appropriate mathematical formulae with words just not to offend the illiterate
  * try to avoid [proprietary](https://en.wikipedia.org/wiki/Proprietary_software) programs and favour [free and open source](https://en.wikipedia.org/wiki/Free_and_open-source_software) ones.
 
 dnl You can ask for help in our mailing list at <wasora@seamplex.com>. There is a community of engineers willing to help you in case you get in trouble with the repositories, the script or the input files.
 
+
 About your favourite FEM program, ask yourself these two questions:
 
  1. Does your favourite FEM program’s manual say what the program does?
@@ -1137,5 +1160,32 @@ About your favourite FEM program, ask yourself these two questions:
 
 And finally, make sure that at the end of the journey from college theory to an actual engineering problem your conscience is clear knowing that there exists a report with your signature on it. That is why we all went to college in the first place.
 
+### Online stuff {#sec:online}
+
+Here is a list of sub-problems and stuff to play with.
+
+ * The pendulum-swing video from [@sec:intro]
+    - <https://youtu.be/Q-lKK4A2OzA>
+ * “On convergence of linearized stresses in an infinite pipe computed using the finite element method” from [@sec:infinite-pipe;@sec:linearity;@sec:infinite-pipe;@sec:infinite-pipe-fem]
+    - <https://www.seamplex.com/fino/doc/pipe-linearized/> 
+ * The three cubes from [@sec:linearity]
+    - case A: pure normal loads (<https://caeplex.com/p/d8fe>)
+    - case B: pure shear loads (<https://caeplex.com/p/b494>)
+    - case C: the combination of A & B (<https://caeplex.com/p/9899>)
+ * The animations of natural oscillations in [@fig:modes]
+    - <https://www.seamplex.com/docs/nafems4/mode1.webm>
+    - <https://www.seamplex.com/docs/nafems4/mode2.webm>
+    - ...
+ * The animations of natural oscillations in [@fig:case-mode]
+    - <https://www.seamplex.com/docs/nafems4/case-mode1.webm>
+    - <https://www.seamplex.com/docs/nafems4/case-mode2.webm>
+    - ...
+ * The parametric tee repository from [@sec:parametric]
+   - <https://github.com/seamplex/tee>
+ * The environmental fatigue sample problem repository from [@sec:in-air;@sec:in-water]
+   - <https://github.com/seamplex/cufen>
+
+See <https://www.seamplex.com/nafems> for new material, updated links and the full version of this case with many more details about the case and the associated mathematics.
+
 
 # Concluding Remarks

nafems4-eps.md

@@ -9,13 +9,18 @@
 
 # CS4: Fatigue Assessment in Nuclear Piping
 
-## From college theory to an actual engineering problem
+## From college theory to an actual engineering problem {#sec:intro}
 
 First of all, please take this text as a written chat between you an me, i.e. an average engineer that has already taken the journey from college to performing actual engineering works using finite element analysis and has something to say about it. Picture yourself in a coffee bar, talking and discussing concepts and ideas with me. Maybe needing to go to a blackboard (or notepad?). Even using a tablet to illustrate some three-dimensional results. But always as a chat between colleagues.
 
 Please also note that I am not a mechanical engineer, although I shared many undergraduate courses with some of them. I am a nuclear engineer with a strong background in mathematics and computer programming. I went to college between 2002 and 2008. Probably a lot of things have changed since then---at least that is what these “millennial” guys and girls seem to be boasting about---but chances are we all studied solid mechanics and heat transfer with a teacher using a piece of chalk on a blackboard while we as students wrote down notes with pencils on paper sheets. And there is really not much that one can do with pencil and paper regarding mechanical analysis. Any actual case worth the time of an engineer needs to be more complex than an ideal canonical case with analytical solution.
 
-Whether you are a student or a seasoned engineer with many years of experience, you might recall from first year physics courses the introduction of the [simple pendulum](https://en.wikipedia.org/wiki/Pendulum) as a case study\ ([@fig:simple]). You learned that the period does not depend on the hanging mass because the weight and the inertia exactly cancelled each other. Also, that [Galileo](https://en.wikipedia.org/wiki/Galileo_Galilei) said (and [Newton proved](https://www.seamplex.com/wasora/realbook/real-012-mechanics.html)) that for small oscillations the period does not even depend on the amplitude. Someone showed you why it worked this way: because if\ $\sin \theta \approx \theta$ then the motion equations converge to an [harmonic oscillator](https://en.wikipedia.org/wiki/Harmonic_oscillator). It might have been a difficult subject for you back in those days when you were learning physics and calculus at the same time. You might have later studied the [Lagrangian](https://en.wikipedia.org/wiki/Lagrangian_mechanics) and even the [Hamiltonian](https://en.wikipedia.org/wiki/Hamiltonian_mechanics) formulations, added a [parametric excitation](https://en.wikipedia.org/wiki/Parametric_oscillator) and analysed the [chaotic double pendulum](https://www.seamplex.com/wasora/realbook/real-017-double-pendulum.html). But it was probably after college, say when you took your first son to a swing on a windy day ([@fig:hamaca]), that you were faced with a real pendulum worth your full attention. Yes, this is my personal story but it could have easily been yours as well. My point is that the very same distance between what I imagined studying a pendulum was and what I saw that day at the swing (namely that the period _does_ depend on the hanging mass) is the same distance between the mechanical problems studied in college and the actual cases encountered during a professional engineer’s lifetime ([@fig:pipes]). In this regard, I am referring only to technical issues. The part of dealing with clients, colleagues, bosses, etc. which is definitely not taught in engineering schools (you can get a heads up in business schools, but again it would be a theoretical pendulum) is way beyond the scope of both this article and my own capacities.
+Whether you are a student or a seasoned engineer with many years of experience, you might recall from first year physics courses the introduction of the [simple pendulum](https://en.wikipedia.org/wiki/Pendulum) as a case study\ ([@fig:simple]). You learned that the period does not depend on the hanging mass because the weight and the inertia exactly cancelled each other. Also, that [Galileo](https://en.wikipedia.org/wiki/Galileo_Galilei) said (and [Newton proved](https://seamplex.com/wasora/doc/realbook/012-mechanics/)) that for small oscillations the period does not even depend on the amplitude. Someone showed you why it worked this way: because if\ $\sin \theta \approx \theta$ then the motion equations converge to an [harmonic oscillator](https://en.wikipedia.org/wiki/Harmonic_oscillator). It might have been a difficult subject for you back in those days when you were learning physics and calculus at the same time.
+divert(-1)
+You might have later studied the [Lagrangian](https://en.wikipedia.org/wiki/Lagrangian_mechanics) and even the [Hamiltonian](https://en.wikipedia.org/wiki/Hamiltonian_mechanics) formulations, added a [parametric excitation](https://en.wikipedia.org/wiki/Parametric_oscillator) and analysed the [chaotic double pendulum](https://seamplex.com/wasora/doc/realbook/017-double-pendulum). divert(0)
+But it was probably after college, say when you took your first son to a swing on a windy day ([@fig:hamaca]), that you were faced with a real pendulum worth your full attention.
+dnl Yes, this is my personal story but it could have easily been yours as well.
+The very same distance between what I imagined studying a pendulum was and what I saw that day at the swing (namely that the period _does_ depend on the hanging mass) is the same distance between the mechanical problems studied in college and the actual cases encountered during a professional engineer’s lifetime ([@fig:pipes]). In this regard, I am referring only to technical issues. The part of dealing with clients, colleagues, bosses, etc. which is definitely not taught in engineering schools (you can get a heads up in business schools, but again it would be a theoretical pendulum) is way beyond the scope of both this article and my own capacities.
 
 ::::: {#fig:pendulum}
 ![Simple pendulum](simple.eps){#fig:simple width=35%}\  
@@ -28,17 +33,18 @@ A [simple pendulum from college physics courses](https://www.seamplex.com/wasora
 ![College pipe](infinite-pipe.eps){#fig:infinite-pipe width=40%}
 ![Real-life pipe](isometric.eps){#fig:isometric width=58%}
 
-An infinitely-long pressurised thick pipe as taught in college and an isometric drawing of a section of a real-life piping system
+An infinitely-long pressurised thick pipe as taught in college and an isometric drawing of a section of a real-life piping system. The details of the isometric are not important individually but they are included to emphasize the fact that a real engineering drawing is far more complex and has far more information than what a professor can draw in a class blackboard.
 :::::
 
+divert(-1)
 Like the pendulums above, we will be swinging back and forth between a case study about fatigue assessment of piping systems in a nuclear power plant and more generic topics related to finite elements and computational mechanics. These latter regressions will not remain just as abstract theoretical ideas. Not only will they be directly applicable to the development of the main case, but they will also apply to a great deal of other engineering problems tackled with the finite element method (and its cousins, [@sec:formulations]).
 
 \medskip
 
 Finite elements are like magic to me. I mean, I can follow the whole derivation of the equations, from the strong, weak and variational formulations of the equilibrium equations for the mechanical problem (or the energy conservation for heat transfer) down to the [algebraic multigrid](https://en.wikipedia.org/wiki/Multigrid_method) [preconditioner](https://en.wikipedia.org/wiki/Preconditioner) for the inversion of the [stiffness matrix](https://en.wikipedia.org/wiki/Stiffness_matrix) passing through [Sobolev spaces](https://en.wikipedia.org/wiki/Sobolev_space) and the [grid generation](https://en.wikipedia.org/wiki/Mesh_generation). Then I can sit down and program all these steps into a computer, including the [shape functions](https://en.wikipedia.org/wiki/Finite_element_method_in_structural_mechanics#Interpolation_or_shape_functions) and its derivatives, the assembly of the discretised stiffness matrix ([@sec:building]), the numerical solution of the system of equations ([@sec:solving]) and the computation of the gradient of the solution (@sec:stress-computation). Yet, the fact that all these a-priori unconnected steps give rise to pretty pictures that resemble reality is still astonishing to me.
+divert(0)
 
-
-Again, take all this information as coming from a fellow that has already taken such a journey from college’s pencil and paper to real engineering cases involving complex numerical calculations. And developing, in the meantime, both an actual working finite-element [back-end](https://www.seamplex.com/fino) and [front-end](https://www.caeplex.com) from scratch.^[The dean of the engineering school I attended used to say “It is not the same to read than to write manuals, and we should aim at writing them.”]
+Again, take all this information as coming from a fellow that has already taken such a journey from college’s pencil and paper to real engineering cases involving complex numerical calculations. And developing, in the meantime, both an actual working finite-element [back-end](https://www.seamplex.com/fino) and [front-end](https://www.caeplex.com) from scratch.dnl^[The dean of the engineering school I attended used to say “It is not the same to read than to write manuals, and we should aim at writing them.”]
 
 
 ### Tips and tricks
@@ -46,8 +52,10 @@ Again, take all this information as coming from a fellow that has already taken
 There are some useful tricks that come handy when trying to solve a mechanical problem. Throughout this text, I will try to tell you some of them.
 
 
-One of the most important ones to use your _imagination_. You will need a lot of imagination to “see” what it is actually going on when analysing an engineering problem. This skill comes from my background in nuclear engineering where I had no choice but to imagine a [positron-electron annihilation](https://en.wikipedia.org/wiki/Electron%E2%80%93positron_annihilation) or an [spontaneous fission](https://en.wikipedia.org/wiki/Spontaneous_fission). But in mechanical engineering, it is likewise important to be able to imagine how the loads “press” one element with the other, how the material reacts depending on its properties, how the nodal displacements generate stresses (both normal and shear), how results converge, etc. And what these results actually mean besides the pretty-coloured figures.^[A former boss once told me “I need the CFD” when I handed in some results. I replied that I did not do computational fluid-dynamics but computed the neutron flux kinetics within a nuclear reactor core. He joked “I know, what I need are the _Colors For Directors_, those pretty-coloured figures along with your actual results to convince the managers.”]
-This journey will definitely need your imagination. We will take a look at equations, numbers, plots, schematics, CAD geometries, 3D\ views, etc. Still, when the theory says “thermal expansion produces normal stresses” you have to picture in your head three little arrows pulling away from the same point in three directions, or whatever mental picture you have about what you understand thermally-induced stresses are. What comes to your mind when someone says that out of the nine elements of the [stress tensor](https://en.wikipedia.org/wiki/Cauchy_stress_tensor) ([@sec:tensor]) there are only six that are independent? Whatever it is, try to practice that kind of graphical thoughts with every new concept. Nevertheless, there will be particular locations of the text where imagination will be most useful. I will bring the subject up then and again.
+One of the most important ones to use your _imagination_. You will need a lot of imagination to “see” what it is actually going on when analysing an engineering problem. This skill comes from my background in nuclear engineering where I had no choice but to imagine a [positron-electron annihilation](https://en.wikipedia.org/wiki/Electron%E2%80%93positron_annihilation) or an [spontaneous fission](https://en.wikipedia.org/wiki/Spontaneous_fission). But in mechanical engineering, it is likewise important to be able to imagine how the loads “press” one element with the other, how the material reacts depending on its properties, how the nodal displacements generate stresses (both normal and shear), how results converge, etc.
+dnl And what these results actually mean besides the pretty-coloured figures.^[A former boss once told me “I need the CFD” when I handed in some results. I replied that I did not do computational fluid-dynamics but computed the neutron flux kinetics within a nuclear reactor core. He joked “I know, what I need are the _Colors For Directors_, those pretty-coloured figures along with your actual results to convince the managers.”]
+
+This journey will definitely need your imagination. We will peek a little bit into equations, numbers, plots, schematics, CAD geometries, 3D\ views, etc. Still, when the theory says “thermal expansion produces normal stresses” you have to picture in your head three little arrows pulling away from the same point in three directions, or whatever mental picture you have about what you understand thermally-induced stresses are. What comes to your mind when someone says that out of the nine elements of the [stress tensor](https://en.wikipedia.org/wiki/Cauchy_stress_tensor) ([@sec:tensor]) there are only six that are independent? Whatever it is, try to practice that kind of graphical thoughts with every new concept. Nevertheless, there will be particular locations of the text where imagination will be most useful. I will bring the subject up then and again.
 
 
 Another heads up is that we will be digging into some mathematics. Probably they would be simple and you would deal with them very easily. But chances are you do not like equations. No problem! Just ignore them for now. Read the text skipping them, it should work as well.
@@ -67,13 +75,13 @@ After further years passed by, engineers (probably the same people that forked s
 
 ![Three-dimensional CAD model for the piping system of\ [@fig:isometric]](real-piping.png){#fig:real-life width=75%}
 
+
 ### Nuclear reactors
 
 In each of the countries that have at least one nuclear power plant there exists a national regulatory body who is responsible for allowing the owner to operate the reactor. These operating licenses are time-limited, with a range that can vary from 25 to 60 years, depending on the design and technology of the reactor. Once expired, the owner might be entitled to an extension, which the regulatory authority can accept provided it can be shown that a certain (and very detailed) set of safety criteria are met. One particular example of requirements is that of fatigue in pipes, especially those that belong to systems that are directly related to the reactor safety.
 
 ### Pressurised pipes
 
-
 How come that pipes are subject to fatigue? Well, on the one hand and without getting into many technical details, the most common nuclear reactor design uses liquid water as coolant and moderator. On the other hand, nuclear power plants cannot by-pass the thermodynamics of the [Carnot cycle](https://en.wikipedia.org/wiki/Carnot_cycle), and in order to maximise the efficiency of the conversion of the energy stored in the uranium nuclei into electricity we need to reach temperatures as high as possible. So, if we want to have liquid water in the core as hot as possible, we need to increase the pressure. The limiting temperature and pressure are given by the [critical point of water](https://en.wikipedia.org/wiki/Critical_point_(thermodynamics)), which is around 374ºC and 22\ MPa. It is therefore expected to have temperatures and pressures near those values in many systems of the plant, especially in the primary circuit and those that directly interact with it, such as pressure and inventory control system, decay power removal system, feedwater supply system, emergency core-cooling system, etc.
 
 Nuclear power plants are not always working at 100% of their maximum power capacity. They need to be maintained and refuelled, they may undergo operational (and some incidental) transients, they might operate at a lower power due to load following conditions, etc. These transient cases involve changes both in temperatures and in pressures that the pipes are subject to, which in turn give rise to changes in the stresses within the pipes. As the transients are postulated to occur cyclically during a number of times during the life-time of the plant (plus its extension period), mechanical fatigue in these piping systems may arise, especially at the interfaces between materials with different thermal expansion coefficients.
@@ -138,7 +146,7 @@ What does this all have to do with mechanical engineering? Well, once we know wh
 
 ### An infinitely-long pressurised pipe {#sec:infinite-pipe}
 
-Let us proceed to our second step, and consider the infinite pipe subject to uniform internal pressure already introduced in\ [@fig:infinite-pipe]. Actually, we are going to solve the mechanical problem on an infinite hollow cylinder, which looks like pipe. This case is usually tackled in college courses, and chances are you already solved it. In fact, the first (and simpler) problem is the “thin cylinder problem.” Then, the “thick cylinder problem” is introduced (the one we solve below), which is slightly more complex. Nevertheless, it has an [analytical solution](https://www.seamplex.com/fino/doc/pipe-linearized/). For the present case, let us consider an infinite pipe (i.e. a hollow cylinder) of internal radius $a$ and external radius $b$ with uniform mechanical properties---Young’s modulus $E$ and Poisson’s ratio $\nu$---subject to an internal uniform pressure $p$.
+Let us proceed to our second step, and consider the infinite pipe subject to uniform internal pressure already introduced in\ [@fig:infinite-pipe]. Actually, we are going to solve the mechanical problem on an infinite hollow cylinder, which looks like pipe. This case is usually tackled in college courses, and chances are you already solved it. In fact, the first (and simpler) problem is the “thin cylinder problem.” Then, the “thick cylinder problem” is introduced (the one we solve below), which is slightly more complex. Nevertheless, it has an analytical solution.^[A detailed analysis of the analytical solution and how results obtained with finite elements converge with respect to the mesh size can be found in the report “On convergence of linearized stresses in an infinite pipe computed using the finite element method.” See [@sec:online].] For the present case, let us consider an infinite pipe (i.e. a hollow cylinder) of internal radius $a$ and external radius $b$ with uniform mechanical properties---Young’s modulus $E$ and Poisson’s ratio $\nu$---subject to an internal uniform pressure $p$.
 
 What follows is more or less what we are taught in school: some equations with a brief explanation of the results. And then we move on to the next subject.
 
@@ -160,7 +168,7 @@ What does this mean? Well, that overall the whole pipe expands a little bit radi
  1. linearly with the pressure, i.e. twice the pressure means twice the displacement, and
  2. inversely proportional to the Young’s Modulus\ $E$ divided by $1+\nu$, i.e. the more resistant the material, the less radial displacements.
 
-That is how an infinite pipe withstands internal pressure.
+That is how an infinite pipe withstands internal pressure. And that is what we are taught in college, which is actually true by the way!
 
 #### Stresses
 
@@ -333,7 +341,7 @@ Remember the main issue of the fatigue analysis in these systems is to analyse w
 
 
 ::::: {#fig:weldolet}
-![Overall view](weldolet-cad1.png){#fig:weldolet-cad width=65%}
+![Overall view](weldolet-cad1-commented.png){#fig:weldolet-cad width=65%}
 
 ![Unstructured tetrahedra-based grid](weldolet-mesh2.png){#fig:weldolet-mesh width=85%}
 
@@ -344,7 +352,7 @@ CAD model of a piping system with a 3/4-inch weldolet-type fork (stainless steel
 ![Location of six SCLs defined to analyse fatigue around a junction.](weldolet-scls-n.png){#fig:weldolet-scls width=75%}
 
 
-On the one hand, a reasonable number of nodes (it is the number of nodes that defines the problem size, not the number of elements as discussed in [@sec:elements-nodes]) in order to get a decent grid is around 200k for each system. On the other hand, solving a couple of dozens of transient heat transfer problems (which we cannot avoid due to the large thermal inertia of the pipes) during a few thousands of seconds over a couple hundred of thousands of nodes might take more time and storage space to hold the results than we might expect.
+On the one hand, a reasonable number of nodes (it is the number of nodes that defines the problem size, not the number of elements) in order to get a decent grid is around 200k for each system. On the other hand, solving a couple of dozens of transient heat transfer problems (which we cannot avoid due to the large thermal inertia of the pipes) during a few thousands of seconds over a couple hundred of thousands of nodes might take more time and storage space to hold the results than we might expect.
 
 There is a wonderful essay by [Isaac Asimov](https://en.wikipedia.org/wiki/Isaac_Asimov) called [“The Relativity of Wrong”](https://en.wikipedia.org/wiki/The_Relativity_of_Wrong) where he introduces the idea that even if something cannot be computed exactly, there are different levels of error. For instance, believing that the Earth is a sphere is less wrong than believing that the Earth is flat, but wrong nonetheless, since it really deviates from a perfect sphere and resembles more an oblate spheroid.
 
@@ -358,7 +366,7 @@ We can then merge this idea by Asimov with an adapted version of the [Saint-Vena
 
 As an example, let us consider the system depicted in\ [@fig:valve] where there is a stainless-carbon steel interface at the discharge of the valves. Instead of solving the transient heat-conduction problem with the internal temperature of the pipes equal to the temperature of the water in the reference transient condition of the power plant and an external condition of natural convection to the ambient temperature in the whole mesh, a reduced model consisting of half of one of the two valves and a small length of the pipes at both the valve inlet and outlet is used. Once the temperature distribution\ $\hat{T}(\mathbf{x},t)$ for each time is obtained in the reduced mesh ([@fig:valve-temp], which has the origin at the centre of the valve), the actual temperature distribution\ $T(\mathbf{x},t)$ is computed by an algebraic generalisation of $\hat{T}(\mathbf{x},t)$ in the full coordinate system. As stated above, those locations which are not covered by the reduced model are generalised with a time-dependent uniform temperature which is the average of the inner and outer temperatures at the inlet and outlet of the reduced mesh.
 
-![Reduced mesh around a valve refined around the interface where the transient heat conduction problem is solved.](valve-temp.png){#fig:valve-temp width=80%}
+![Mesh of half-a-symmetric-valve refined around the interface where the transient heat conduction problem is solved.](valve-temp-commented.png){#fig:valve-temp width=80%}
 
 
 Note that there is no need to have a one-to-one correspondence between the elements from the reduced mesh with the elements from the original one. Actually, the reduced mesh contains first-order elements whilst the former has second-order elements. Also the grid density is different. Nevertheless, the finite-element solver [Fino](https://www.seamplex.com/fino) used to solve both the heat and the mechanical problems, allows to read functions of space and time defined over one mesh and continuously evaluate and use them into another one even if the two grids have different elements, orders or even dimensions. In effect, in the system from\ [@fig:real-life] the material interface is between a orifice plate made in stainless steel that is welded to a carbon-steel pipe ([@fig:real]). The thermal problem can be modelled using a two-dimensional axi-symmetric grid and then generalised to the full three-dimensional mesh using the algebraic manipulation capabilities provided by [Fino](https://www.seamplex.com/fino) (actually by [wasora](https://www.seamplex.com/wasora)) as shown in\ [@fig:real].
@@ -366,16 +374,17 @@ Note that there is no need to have a one-to-one correspondence between the eleme
 dnl ![Full 3D mesh, carbon steel is magenta, stainless steel is green.](real-mesh2.png){#fig:real-mesh2 width=70%}
 dnl ![Temperature distribution for a certain time within the transient computed on a reduced two-dimensional axi-symmetric mesh modelling half the orifice plate and a length of the carbon pipe.](real-temp.png){#fig:real-temp width=34%}
 
-![Temperature distribution for a certain time within the transient computed on a reduced two-dimensional axi-symmetric mesh modelling half the orifice plate and a length of the carbon pipe and Generalisation of the temperature to the full three-dimensional grid.](real-gen.png){#fig:real width=85%}
+![Temperature distribution for a certain time within the transient computed on a reduced two-dimensional axi-symmetric mesh modelling half the orifice plate and a length of the carbon pipe and generalisation of the temperature to the full three-dimensional grid.](real-gen.png){#fig:real width=85%}
 
 
 ### Seismic loads {#sec:seismic}
 
-Before considering the actual mechanical problem that will give us the stress tensor at the SCLs, and besides needing to solve the transient thermal problem to get the temperature distributions, we need to address the loads that arise due to a postulated earthquake during a certain part of the operational transients. The full computation of a mechanical transient problem using the earthquake time-dependent displacements is off the table for two reasons. First, because again the computation would take more time than we might have to deliver the report. And secondly and more importantly, because civil engineers do not compute earthquakes in the time domain but in the frequency domain using the [response spectrum method](https://en.wikipedia.org/wiki/Response_spectrum), provided the case is linear (see\ [@sec:linearity]). Time to revisit our [Fourier transform](https://en.wikipedia.org/wiki/Fourier_transform) exercises from undergraduate math courses.
 
+Every nuclear power plant is designed to withstand earthquakes. Of course, not all plants need the same level of reinforcements. Those built in large quiet plains will be, seismically speaking, cheaper than those located in geologically active zones. Keep in mind that all the 54 Japanese nuclear power plants did structurally resist the 2011 earthquake, and all of the reactors were safely shut down. What actually happened in [Fukushima](http://www.world-nuclear.org/information-library/safety-and-security/safety-of-plants/fukushima-accident.aspx) is that one hour after the main shake, a 14-metre tsunami splashed on the coast, jumping over the 9-metre defences and flooding the emergency [Diesel generators](https://en.wikipedia.org/wiki/Diesel_generator) that provided power to the pumps in charge of removing the remaining [decay power](https://en.wikipedia.org/wiki/Decay_heat) from the already-stopped [reactor core](https://en.wikipedia.org/wiki/Nuclear_reactor_core).
+
+divert(-1)
 #### Earthquake spectra
 
-Every nuclear power plant is designed to withstand earthquakes. Of course, not all plants need the same level of reinforcements. Those built in large quiet plains will be, seismically speaking, cheaper than those located in geologically active zones. Keep in mind that all the 54 Japanese nuclear power plants did structurally resist the 2011 earthquake, and all of the reactors were safely shut down. What actually happened in [Fukushima](http://www.world-nuclear.org/information-library/safety-and-security/safety-of-plants/fukushima-accident.aspx) is that one hour after the main shake, a 14-metre tsunami splashed on the coast, jumping over the 9-metre defences and flooding the emergency [Diesel generators](https://en.wikipedia.org/wiki/Diesel_generator) that provided power to the pumps in charge of removing the remaining [decay power](https://en.wikipedia.org/wiki/Decay_heat) from the already-stopped [reactor core](https://en.wikipedia.org/wiki/Nuclear_reactor_core).
 
 Back to our case study, the point is that each site where nuclear power plants are built must have a geological study where a postulated design-basis earthquake is to be defined. In other words, a theoretical earthquake which the plant ought to withstand needs to be specified. How? By giving a set of three spectra (one for each coordinate direction) giving acceleration as a function of the frequency for each level of the building. That is to say, once the earthquake hits the power plant, depending on soil-structure interactions the energy will shake the building foundations in a way that depends on the characteristics of the earthquake, the soil and the concrete structure. Afterwards, the way the oscillations travel upward and shake each of the mechanical components erected in each floor level depends on the design of the civil structure in a way which is fully determined by the floor response spectra like the ones depicted in\ [@fig:spectrum].
 
@@ -399,8 +408,13 @@ A real continuous solid has infinite modes of oscillation. A discretised one (us
 
 ![Effective accumulated mass fraction up to the $i$-th mode of oscillation.](acumulada.eps){#fig:acumulada}
 
+
+These first modes, shown in\ [@fig:modes], that take up most of the energy are then used to take into account the earthquake load. There are several ways of performing this computation, but the ASME\ III code states that the method known as SRSS (for Square Root of Sum of Squares) can be used. This method mixes the eigenvectors with the floor response spectra through the eigenvalues and gives a spatial (actually nodal) distribution of three accelerations (one for each direction) that, when multiplied by the density of the material, give a vector of a distributed force (in units of Newton per cubic millimetre for example) which is statically equivalent to the load coming from the postulated earthquake.
+divert(0)
+
+
 ::::: {#fig:modes}
-![$f_1 \approx 30$\ [Hz](https://www.seamplex.com/docs/nafems4/mode1.webm)](mode1.png){width=48%}\  
+![$f_1 \approx 30$\ [Hz](https://www.seamplex.com/docs/nafems4/mode1.webm)](mode1-commented.png){width=48%}\  
 ![$f_2 \approx 60$\ [Hz](https://www.seamplex.com/docs/nafems4/mode2.webm)](mode2.png){width=48%}
 
 ![$f_3 \approx 70$\ [Hz](https://www.seamplex.com/docs/nafems4/mode3.webm)](mode3.png){width=48%}\  
@@ -409,12 +423,9 @@ A real continuous solid has infinite modes of oscillation. A discretised one (us
 ![$f_5 \approx 100$\ [Hz](https://www.seamplex.com/docs/nafems4/mode5.webm)](mode5.png){width=48%}\  
 ![$f_6 \approx 130$\ [Hz](https://www.seamplex.com/docs/nafems4/mode6.webm)](mode6.png){width=48%}
 
-First six natural oscillation modes for a piping section between axial supports. Blue-red color scale shows magnitude of deformation, as compared to the original geometry (grey). Animations can be seen at links like <https://www.seamplex.com/docs/nafems4/mode1.webm>.
+First six natural oscillation modes for a piping section between axial supports. Blue-red colour scale shows magnitude of deformation, as compared to the original geometry (grey). Full animations available online ([@sec:online]).
 :::::
 
-
-These first modes, shown in\ [@fig:modes], that take up most of the energy are then used to take into account the earthquake load. There are several ways of performing this computation, but the ASME\ III code states that the method known as SRSS (for Square Root of Sum of Squares) can be used. This method mixes the eigenvectors with the floor response spectra through the eigenvalues and gives a spatial (actually nodal) distribution of three accelerations (one for each direction) that, when multiplied by the density of the material, give a vector of a distributed force (in units of Newton per cubic millimetre for example) which is statically equivalent to the load coming from the postulated earthquake.
-
 ::::: {#fig:acceleration}
 ![$a_x$](ax.png){width=80%}
 
@@ -422,10 +433,11 @@ These first modes, shown in\ [@fig:modes], that take up most of the energy are t
 
 ![$a_z$](az.png){width=80%}
 
-The equivalent accelerations for the piping section of [@fig:modes] for the spectra of\ [@fig:spectrum]
+Equivalent accelerations for a certain piping section
 :::::
 
-The ASME code says that these accelerations (depicted in [@fig:acceleration]) ought to be applied twice: once with the original sign and once with all the elements with the opposite sign, as SRSS “looses signs.” The application of each of these equivalent loads should last two seconds in the original time domain.
+Since the computation of the loads that a certain earthquake gives rise to would add a significant amount of complexity to this already-complex case study, the details are skipped. In any case, we need to compute the first natural modes of oscillation of the piping section under study ([@fig:modes]) and then, after some cumbersome maths, obtain an statically-equivalent load distribution. In effect, [@fig:acceleration] shows a sample distribution of acceleration distribution within a certain piping system which, when multiplied by the density, give a load distribution which is statically equivalent to the dynamic solicitations created by the earthquake.
+The ASME code says that these accelerations ought to be applied twice: once with the original sign and once with all the elements with the opposite sign, as SRSS “looses signs.” The application of each of these equivalent loads should last two seconds in the original time domain.
 
 
 ### Linearity (not yet linearisation) {#sec:linearity}
@@ -436,17 +448,17 @@ Let us jump out of our nuclear piping problem and step back into the general fin
 
 \medskip
 
-Time for both of us to make an experiment. Grab your favourite FEM program for the first time (remember mine is [CAEplex](https://caeplex.com)) and create a 1mm $\times$ 1mm $\times$ 1mm cube. Set any values for the Young’s Modulus and Poisson ratio as you want. I chose\ $E=200$\ GPa and\ $\nu=0.28$. Restrict the three faces pointing to the negative axes to their planes, i.e.
+Time for both of us to make an experiment. Grab your favourite FEM program for the first time (remember mine is [CAEplex](https://caeplex.com), which can also be accessed through [Onshape](https://www.onshape.com/cad-blog/partner-spotlight-using-caeplex-for-finite-element-analysis-in-onshape)) and create a 1mm $\times$ 1mm $\times$ 1mm cube. Set any values for the Young’s Modulus and Poisson ratio as you want. I chose\ $E=200$\ GPa and\ $\nu=0.28$. Restrict the three faces pointing to the negative axes to their planes, i.e.
 
  * in face “left” ($x<0$), set null displacement in the $x$ direction ($u=0$),
  * in face “front” ($y<0$), set null displacement in the $y$ direction ($v=0$),
  * in face “bottom” ($z<0$), set null displacement in the $z$ direction ($w=0$).
 
-Now we are going to create and compare three load cases:
+Now we are going to create and compare three load cases:^[The provided links lead to FEA cases which are fully usable directly from the web browser, i.e. “finite elements in the cloud.”]
 
- a. Pure normal loads (<https://caeplex.com/p/d8f>)
+ a. Pure normal loads (<https://caeplex.com/p/d8fe>)
  b. Pure shear loads (<https://caeplex.com/p/b494>)
- c. The combination of A & B (<https://caeplex.com/p/989>)
+ c. The combination of A & B (<https://caeplex.com/p/9899>)
 
 The loads in each cases are applied to the three remaining faces, namely “right” ($x>0$), “back” ($y>0$) and “top,” ($z>0$). Their magnitude in Newtons are:
 
@@ -477,8 +489,9 @@ $$
 
 
 ::::: {#fig:cube}
-![Case B, pure-shear loads](cube-shear.png){#fig:cube-shear width=48%}\  
-![Case C, normal plus shear loads](cube-full.png){#fig:cube-full width=48%}
+![Case B, pure-shear loads](cube-shear.png){#fig:cube-shear width=65%}
+
+![Case C, normal plus shear loads](cube-full.png){#fig:cube-full width=65%}
 
 Spatial distribution of principal stress\ 3 for cases\ B and\ C. If linearity applied, case\ C would be equal to case\ B plus a constant
 :::::
@@ -623,14 +636,14 @@ The first thing that has to be said is that, as with any interesting problem, th
  10. etc.
 
 ::::: {#fig:quarter}
-![Structured second-order incomplete hexahedra](quarter-struct.png){#fig:cube-struct width=48%}\  
-![Unstructured second-order tetrahedra](quarter-caeplex.png){#fig:quarter-caeplex width=48%}
+![Structured second-order hexahedra](quarter-struct.png){#fig:cube-struct width=50%}\  
+![Unstructured second-order tetrahedra](quarter-caeplex.png){#fig:quarter-caeplex width=50%}
 
 Two of the hundreds of different ways the infinite pressurised pipe can be solved using FEM. The axial displacement at the ends is set to zero, leading to a [plane-strain](https://en.wikipedia.org/wiki/Plane_stress#Plane_strain_(strain_matrix)) condition
 :::::
 
 
-You can get both the exponential nature of each added bullet and how easily we can add new further choices to solve a FEM problem. And each of these choices will reveal you something about the nature of either the mechanical problem or the numerical solution. It is not possible to teach any possible lesson from every outcome in college, so you will have to learn them by yourself getting your hands at them. I have already tried to address the particular case of the infinite pipe in a [recent report](https://www.seamplex.com/fino/doc/pipe-linearized/)^[<https://www.seamplex.com/fino/doc/pipe-linearized/>] that is worth reading before carrying on with this article.
+You can get both the exponential nature of each added bullet and how easily we can add new further choices to solve a FEM problem. And each of these choices will reveal you something about the nature of either the mechanical problem or the numerical solution. It is not possible to teach any possible lesson from every outcome in college, so you will have to learn them by yourself getting your hands at them. I have already tried to address the particular case of the infinite pipe in a [recent report](https://www.seamplex.com/fino/doc/pipe-linearized/)^[The aforementioned “On convergence of linearized stresses in an infinite pipe computed using the finite element method.”] that is worth reading before carrying on with this article.
 The main conclusions of the report are:
 
 dnl * Engineering problems ought not to be solved using black-boxes (i.e. privative software whose source code is not freely available)---more on the subject below in\ [@sec:two-materials].
@@ -640,12 +653,12 @@ dnl * There are no shear stresses, so these three stresses are also the principa
 dnl * Analytical expressions for the membrane and membrane plus bending stresses along any radial SCL can be obtained.
 dnl * The spatial domain can be discretised using linear or higher-order elements. In particular first and second-order elements have been used in the report.
  * For the same number of elements, second-order grids need more nodes than linear ones, although they can better represent curved geometries.
- * The discretised problem size depends on the number of nodes and not on the number of elements---more on the subject below in\ [@sec:elements-nodes].
+dnl * The discretised problem size depends on the number of nodes and not on the number of elements---more on the subject below in\ [@sec:elements-nodes].
 dnl * The finite-element results for the displacements are similar to the analytical solution, with second-order grids giving better results for the same number of elements (this is expected as they involved far more nodes).
  * The three stress distributions computed with the finite-element give far more reasonable results for the second-order case than for the first-order grid. This is qualitatively explained by the fact that  first-order tetrahedra have uniform derivatives and such the elements located in both the external and external faces represent the stresses not at the actual faces but at the barycentre of the elements.
  * Membrane stresses converge well for both the first and second-order cases because they represent a zeroth-order moment of the stress distribution and the excess and defect errors committed by the internal and external elements approximately cancel out.
  * Membrane plus bending stresses converge very poorly with linear elements because the excess and defect errors do not cancel out because it is a first-order moment of the stress distribution.
- * The computational effort to solve a given problem, namely the CPU time and the needed RAM ([@fig:error-vs-cpu]) depend non-linearly on various factors, but the most important one is the problem size which is three times the number of nodes in the grid---more on the subject below in\ [@sec:elements-nodes].
+ * The computational effort to solve a given problem, namely the CPU time and the needed RAM ([@fig:error-vs-cpu]) depend non-linearly on various factors, but the most important one is the problem size which is three times the number of nodes in the grid dnl---more on the subject below in\ [@sec:elements-nodes].
  * The error with respect to the analytical solutions as a function of the CPU time needed to compute the membrane stress is similar for both first and second-order grids. But for the computation of the membrane plus bending stress ([@fig:error-MB-vs-cpu]), first-order grids give very poor results compared to second-order grids for the same CPU time.
 
 ::::: {#fig:error-vs-cpu}
@@ -657,6 +670,7 @@ Error in the computation of the linearised stresses vs. CPU time needed to solve
 
 An additional note should be added. The FEM solution, which not only gives the nodal displacements but also a method to interpolate these values inside the elements, does not fully satisfy the original equilibrium equations at every point (i.e. the strong formulation). It is an approximation to the solution of the [weak formulation](https://en.wikipedia.org/wiki/Weak_formulation) that is close (measured in the vector space spanned by the [shape functions](https://www.quora.com/What-is-a-shape-function-in-FEM)) to the real solution. Mechanically, this means that the FEM solution leads only to nodal equilibrium but not point-wise equilibrium.
 
+divert(-1)
 ### Elements, nodes and CPU {#sec:elements-nodes}
 
 dnl The last two bullets above lead to an issue that has come many times when discussing the issue of convergence with respect to the mesh size with other colleagues. There apparently exists a common misunderstanding that the number of elements is the main parameter that defines how complex a FEM model is. This is strange, because even in college we are taught that the most important parameter is the _size_ of the stiffness matrix, which is three times (for 3D problems with the [displacement-based formulation](http://web.mit.edu/kjb/www/Books/FEP_2nd_Edition_4th_Printing.pdf) formulation) the number of _nodes_.
@@ -732,7 +746,7 @@ Now proceed to picturing the general three-dimensional cases with unstructured t
 Detailed mathematics show that the location where the derivatives of the interpolated displacements are closer to the real (i.e. the analytical ones in problems that have it) solution are the elements’ [Gauss points](https://en.wikipedia.org/wiki/Gaussian_quadrature). Even better, the material properties at these points are continuous (they are usually uniform but they can depend on temperature for example) because, unless we are using weird elements, there are no material interfaces inside elements. But how to manage a set of stresses given at the Gauss points instead of at the nodes? Should we use one mesh for the input and another one for the output? What happens when we need to know the stresses on a surface and not just in the bulk of the solid? There are still no one-size-fits-all answers. There is a very interesting [blog post](https://nickjstevens.netlify.com/post/2019/stress-singularities/) by Nick Stevens that addresses the issue of stresses computed at sharp corners. What does your favourite FEM program do with such a case?
 
 In any case, this step takes a non-negligible amount of time. The most-common approach, i.e. the node-averaging method is driven mainly by the number of nodes of course. So all-in-all, these are the reasons to use the number of nodes instead of the numbers of elements as a basic parameter to measure the complexity of a FEM problem.
-
+divert(0)
 
 ## Adding complexity: the truth is out there
 
@@ -861,7 +875,7 @@ Do you now see the added value of training throw-ins with watermelons? We might
  * using distributed forces from earthquakes,
  * etc.
 
-Most of the time at college we would barely do what is needed to be approved in one course and move on to the next one. If you have the time and consider a career related to finite-element analysis, please do not. Clone the repository^[<https://github.com/seamplex/tee>] with the input files for [Fino](https://www.seamplex.com/fino) and start playing. If you are stuck, do not hesitate asking for help in [wasora’s mailing list](https://www.seamplex.com/lists.html).
+Most of the time at college we would barely do what is needed to be approved in one course and move on to the next one. If you have the time and consider a career related to finite-element analysis, please do not. Clone the repository ([@sec:online]) with the input files for [Fino](https://www.seamplex.com/fino) and start playing. If you are stuck, do not hesitate asking for help in [wasora’s mailing list](https://www.seamplex.com/lists.html).
 
 One further detail: it is always a sane check to try to explain the numerical results based on physical reasoning (i.e. “with your fingers”) as we did two paragraphs above. Most of the time you will be solving problems whilst already knowing what the result would (or ought to) be.
 
@@ -872,7 +886,7 @@ A fellow mechanical engineer who went to the same high school I did, who went to
 
 So here we are again with the case study where we have to compute the linearised stresses at certain SCLs located near the interface between two different kinds of steels during operational and incidental transients of the plant. The first part is then to “bake” the pipes, modelling a thermal transient with time-dependent temperature or convection boundary conditions (it depends on the system). This steps gives a time and space-dependent temperature\ $T(x,y,z,t)$.
 
-Then we proceed to “shake” the pipes, i.e. to compute the natural frequencies and associated oscillations modes. Employing the floor response spectra and combining the individual contributions with the SRSS method discussed in section\ [@sec:seismic], we obtain a distributed load vector\ $\mathbf{f}(x,y,z)$ which is statically equivalent to the design earthquake.
+We then proceed to “shake” the pipes. That is to say, we obtain a distributed load vector\ $\mathbf{f}(x,y,z)$ which is statically equivalent to the design earthquake.
 
 Finally we attempt to “break” the pipes successively solving many steady-state elastic problems for different times\ $t$ of the operational transient. Some remarks about this step:
 
@@ -908,7 +922,11 @@ If you notice, the complex dependence of the stiffness matrix\ $K$ can be reduce
 
 $$ K(\mathbf{x}) \cdot \mathbf{u}(\mathbf{x}) = \mathbf{b}(\mathbf{x})$$
 
-And this last equation is linear in\ $\mathbf{u}$. In effect, the discretisation step means to integrate over\ $\mathbf{x}$. As\ $K$, $\mathbf{u}$ and\ $\mathbf{b}$ depend only on\ $\mathbf{x}$, then after integration one gets just numbers with the matrix representation of\ [@eq:kub]. Again, you can either trust me, ask a teacher or go through with the maths (in increasing order of recommendation).
+And this last expression is linear in\ $\mathbf{u}$! In effect, the discretisation step means to integrate over\ $\mathbf{x}$. As\ $K$, $\mathbf{u}$ and\ $\mathbf{b}$ depend only on\ $\mathbf{x}$, then after integration one gets just numbers inside $K$ and $\mathbf{b}$. Again, you can either, in increasing order of recommendation:
+
+ 1. trust me,
+ 2. ask a teacher, or
+ 3. go through with the maths.
 
 
 ::::: {#fig:case}
@@ -943,7 +961,7 @@ Temperature distribution for a certain instant of the transient, computed in the
 ![$f_1 \approx 930$\ [Hz](https://www.seamplex.com/docs/nafems4/case-mode9.webm)](case-mode9.png){#fig:case-mode-9 width=30%}
 
 First nine natural modes of oscillation of the piping system subject to the boundary conditions the supports provide.
-Animations can be seen at links like <https://www.seamplex.com/docs/nafems4/case-mode1.webm>.
+Animations are available online ([@sec:online]).
 :::::
 
 dnl ![Floor response spectrum.](case-spectrum.eps){#fig:case-spectrum width=70%}
@@ -955,15 +973,15 @@ To recapitulate, the figures in this section show some partial non-dimensional r
  3. defining the number and locations of the SCLs ([@fig:case-scls] or [@fig:weldolet-scls;@fig:valve])
  4. computing a heat conduction (bake) transient problem with temperatures as a function of time from the operational transient in a simple domain using temperature-dependent thermal conduction coefficients from ASME\ II\ part D ([@fig:case-temp] or [@fig:valve-temp])
  5. generalising the temperature distribution as a function of time to the general domain ([@fig:case-temp2] or [@fig:real])
- 6. performing a modal analysis (shake) on the main domain to obtain the main oscillation frequencies and modes ([@fig:case-mode] or [@fig:modes])
- 7. using the floor response spectra ([@fig:spectrum]) and the SRSS method to obtain a distributed force statically-equivalent to the earthquake load ([@fig:case-acceleration] or [@fig:acceleration])
+ 6. performing a modal analysis (shake) on the main domain to obtain the main oscillation frequencies and modes ([@fig:case-mode] dnl or [@fig:modes])
+ 7. obtaining a distributed force statically-equivalent to the earthquake load (shake, [@fig:case-acceleration] or [@fig:acceleration])
  8. successively solving the linear elastic problem for different times using the generalised temperature distribution taking into account
      a. the dependence of the Young’s Modulus\ $E$ and the thermal expansion coefficient\ $\alpha$ with temperature,
      b. the thermal expansion effect itself
      c. the instantaneous pressure exerted in the internal faces of the pipes at the time\ $t$ according to the definition of the operational transient
      d. the restriction of the degrees of freedom of those faces, lines or points that correspond to mechanical supports located both within and at the ends of the CAD model
      e. the earthquake load, which according to ASME should be present only during four seconds of the transient: two seconds with one sign and the other two seconds with the opposite sign. This period should be selected to coincide with the instant of highest mechanical stress (conservative computation)
- 9. computing the linearised stresses (membrane and membrane plus bending) at the SCLs combining them as Principal\ 1,  Principal\ 2, Principal\ 3 and Tresca
+ 9. computing the linearised stresses (membrane and membrane plus bending) at the SCLs combining them as Principal\ 1,  Principal\ 2, Principal\ 3 and Tresca^[A deep technical note should be added for discussing the validity of linearising a stress invariant as the principall stresses individually. It should be enough for the present study to keep in mind that in [@sec:in-air] it is the difference of these linearised principal stresses (i.e. Tresca-like) are used to compute the stress amplitude of the periodic load for fatigue assessment.]
  10. juxtaposing these linearised stresses for each time of the transient and for each transient so as to obtain a single time-history of stresses including all the operational and/or incidental transients under study, which is what stress-based fatigue assessment needs (recall\ [@sec:fatigue] and go on to\ [@sec:usage]).
 
 A pretty nice list of steps, which definitely I would not have been able to tackle when I was in college. Would you?
@@ -1009,7 +1027,8 @@ When\ $\text{CUF} < 1$, then the part under analysis can withstand the proposed
 This cryptic paragraph can be better explained by using a clearer example. To avoid using actual sensitive data from a real power plant, let us use the same test case used by both the [US Nuclear Regulatory Commission](https://en.wikipedia.org/wiki/Nuclear_Regulatory_Commission) (in its report NUREG/CR-6909) and the [Electric Power Institute](https://en.wikipedia.org/wiki/Electric_Power_Research_Institute) (report 1025823) called “EAF (Environmentally-Assisted Fatigue) Sample Problem 2-Rev.\ 2 (10/21/2011)”.
 
 
-![A low-alloy steel vessel nozzle (blue) welded to a stainless steel pipe (grey)](axi-inches-3d.png){#fig:axi-inches-3d width=25%}
+![A low-alloy steel vessel nozzle (blue) welded to a stainless steel pipe (grey)](axi-inches-3d.png){#fig:axi-inches-3d width=35%}
+
 
 
 It consists of a typical vessel nozzle with attached piping as illustrated in\ [@fig:axi-inches-3d]. These components are subject to four transients\ $k=1,2,3,4$ that give rise to linearised stress histories (slightly modified according to NB-3216.2) which are given as individual stress values juxtaposed so as to span a time range of about 36,000 seconds ([@fig:nureg1]). As the time steps is not constant, each stress value has an integer index\ $i$ that uniquely identifies it:
@@ -1021,8 +1040,24 @@ It consists of a typical vessel nozzle with attached piping as illustrated in\ [
 |   3   |  6450--9970    |  960--1595  |         20       |
 |   4   |  9970--35971   |  1596--2215 |        100       |
 
-A design-basis earthquake was assumed to occur briefly during one second (sic) at around\ $t=3470$\ seconds, and it is assigned a number of cycles\ $n_e=5$. The detailed stress history for one of the SCLs including both the principal and lineariased stresses, which are already offset following NB-3216.2 so as to have a maximum stress equal to zero, can be found as an appendix in NRC's NUREG/CR-6909 report, or in the repository with the scripts I prepared for you to play with this problem.^[<https://github.com/seamplex/cufen>]
 
+A design-basis earthquake was assumed to occur briefly during one second (sic) at around\ $t=3470$\ seconds, and it is assigned a number of cycles\ $n_e=5$. The detailed stress history for one of the SCLs including both the principal and lineariased stresses, which are already offset following NB-3216.2 so as to have a maximum stress equal to zero, can be found as an appendix in NRC's NUREG/CR-6909 report, or in the repository with the scripts I prepared for you to play with this problem.
+
+To compute the global usage factor, we first need to find all the combinations of local extrema pairs and then sort them in decreasing order of stress difference. For example, the largest stress amplitude is found between $i=447$ and $i=694$. The second one is 447--699. Then 699--1020, 699--899, etc. For each of these pairs, defined by the indexes\ $i_{1,j}$ and $i_{2,j}$, a partial usage factor\ $U_j$ should computed. The stress amplitude\ $S_{\text{alt},j}$ which should be used to enter into the $S$-$N$ curve is
+
+\pagebreak
+
+
+$$
+S_{\text{alt},j} = \frac{1}{2} \cdot k_{e,j} \cdot \left| MB^\prime_{i_{1,j}} - MB^\prime_{i_{2,j}} \right| \cdot \frac{E_\text{SN}}{E(T_{\text{max}_j})}
+$$
+
+\noindent where $k_e$ is a plastic correction factor for large loads (NB-3228.5), $E_\text{SN}$ is the Young’s Modulus used to create the $S$-$N$ curve (provided in the ASME fatigue curves) and\ $E(T_{\text{max}_j})$ is the material’s Young’s Modulus at the maximum temperature within the\ $j$-th interval.
+
+
+We now need to comply with ASME’s obscure note about the number of cycles to assign a proper value of\ $n_j$. Back to the largest pair 447-694, we see that 447 belongs to transient\ #1 which has assigned 20 cycles and 694 belongs to the earthquake with 5 cycles. Therefore $n_1=5$, and the cycles associated to each index are decreased in five. So $i=694$ disappears and the number of cycles associated to $i=447$ are decreased from 20 to 15. The second largest pair is now 447-699, with 15 (because we just spent 5 in the first pair) and 50 cycles respectively. Now $n_2=15$, point 447 disappears and 699 remains with 35 cycles. The next pair is 699-1020, with number of cycles 35 and 20 so $n_3=20$, point 1020 disappears and 699 remains with 15 cycles. And so on, down to the last pair.
+
+Why all these details? Not because I want to teach you how to perform fatigue evaluations just reading this section without resorting to ASME, fatigue books and even other colleagues. It is to show that even though these computation can be made “by hand” (i.e. using a calculator or, God forbids, a spreadsheet) when having to evaluate a few SCLs within several piping systems it is far (and I mean really far) better to automate all these steps by writing a set of scripts. Not only will the time needed to process the information be reduced, but also the introduction of human errors will be minimised and repeatability of results will be assured---especially if working under a [distributed version control](https://en.wikipedia.org/wiki/Distributed_version_control) system such as [Git](https://en.wikipedia.org/wiki/Git). This is true in general, so here is another tip: learn to write scripts to post-process your FEM results (you will need to use a script-friendly FEM program) and you will gain considerable margins regarding time and quality.
 
 ::::: {#fig:nureg}
 ![Full time range of the juxtaposed transients spanning $\approx$ 36,000 seconds](nureg1.eps){#fig:nureg1 width=70%}
@@ -1034,14 +1069,6 @@ A design-basis earthquake was assumed to occur briefly during one second (sic) a
 Time history of the linearised stress\ $\text{MB}_{31}$ corresponding to the example problem from NRC and EPRI. The indexes\ $i$ of the extrema are shown in green (minimums) and red (maximums)
 :::::
 
-
-To compute the global usage factor, we first need to find all the combinations of local extrema pairs and then sort them in decreasing order of stress difference. For example, the largest stress amplitude is found between $i=447$ and $i=694$. The second one is 447--699. Then 699--1020, 699--899, etc. For each of these pairs, defined by the indexes\ $i_{1,j}$ and $i_{2,j}$, a partial usage factor\ $U_j$ should computed. The stress amplitude\ $S_{\text{alt},j}$ which should be used to enter into the $S$-$N$ curve is
-
-$$
-S_{\text{alt},j} = \frac{1}{2} \cdot k_{e,j} \cdot \left| MB^\prime_{i_{1,j}} - MB^\prime_{i_{2,j}} \right| \cdot \frac{E_\text{SN}}{E(T_{\text{max}_j})}
-$$
-where $k_e$ is a plastic correction factor for large loads (NB-3228.5), $E_\text{SN}$ is the Young’s Modulus used to create the $S$-$N$ curve (provided in the ASME fatigue curves) and\ $E(T_{\text{max}_j})$ is the material’s Young’s Modulus at the maximum temperature within the\ $j$-th interval.
-
 ::::: {#fig:cuf}
 ![Reference from NUREG/CR6909](cuf-nrc.png){#fig:cuf-nrc width=100%}
 
@@ -1050,9 +1077,6 @@ where $k_e$ is a plastic correction factor for large loads (NB-3228.5), $E_\text
 Tables of individual usage factors for the NRC/EPRI “EAF Sample Problem 2-Rev.\ 2 (10/21/2011).” One table is taken from a document issued by almost-a-billion-dollar-year-budget government agency from the most powerful country in the world and the other one is from a third-world engineering startup. Guess which is which.
 :::::
 
-We now need to comply with ASME’s obscure note about the number of cycles to assign a proper value of\ $n_j$. Back to the largest pair 447-694, we see that 447 belongs to transient\ #1 which has assigned 20 cycles and 694 belongs to the earthquake with 5 cycles. Therefore\ $n_1=5$, and the cycles associated to each index are decreased in five. So $i=694$ disappears and the number of cycles associated to $i=447$ are decreased from 20 to 15. The second largest pair is now 447-699, with 15 (because we just spent 5 in the first pair) and 50 cycles respectively. Now $n_2=15$, point\ 447 disappears and 699 remains with 35 cycles. The next pair is 699-1020, with number of cycles 35 and 20 so $n_3=20$, point 1020 disappears and 699 remains with 15 cycles. And so on, down to the last pair.
-
-Why all these details? Not because I want to teach you how to perform fatigue evaluations just reading this section without resorting to ASME, fatigue books and even other colleagues. It is to show that even though these computation can be made “by hand” (i.e. using a calculator or, God forbids, a spreadsheet) when having to evaluate a few SCLs within several piping systems it is far (and I mean really far) better to automate all these steps by writing a set of scripts. Not only will the time needed to process the information be reduced, but also the introduction of human errors will be minimised and repeatability of results will be assured---especially if working under a [distributed version control](https://en.wikipedia.org/wiki/Distributed_version_control) system such as [Git](https://en.wikipedia.org/wiki/Git). This is true in general, so here is another tip: learn to write scripts to post-process your FEM results (you will need to use a script-friendly FEM program) and you will gain considerable margins regarding time and quality.
 
 ### In water (NRC’s extension) {#sec:in-water}
 
@@ -1096,7 +1120,7 @@ Tables of individual environmental correction and usage factors for the NRC/EPRI
 
 ## Conclusions
 
-Back in college, we all learned how to solve engineering problems. And once we graduated, we felt we could solve and fix the world (if you did not graduate yet, you will feel it shortly). But there is a real gap between the equations written in chalk on a blackboard (now probably in the form of beamer slide presentations) and actual real-life engineering problems. This chapter introduces a real case from the nuclear industry and starts by idealising the structure such that it has a known analytical solution that can be found in textbooks.  Additional realism is added in stages allowing the engineer to develop an understanding of the more complex physics and a faith in the veracity of the finite-element results where theoretical solutions are not available. Even more, a brief insight into the world of evaluation of stress-life fatigue using such results further illustrates the complexities of real-life engineering analysis. Here is a list of the tips that arose throughout the text:
+Back in college, we all learned how to solve engineering problems. And once we graduated, we felt we could solve and fix the world (if you did not graduate yet, you will feel it shortly). But there is a real gap between the equations written in chalk on a blackboard (now probably in the form of beamer slide presentations) and actual real-life engineering problems. This chapter introduces a real case from the nuclear industry and starts by idealising the structure such that it has a known analytical solution that can be found in textbooks.  Additional realism is added in stages allowing the engineer to develop an understanding of the more complex physics and a faith in the veracity of the finite-element results where theoretical solutions are not available. Even more, a brief insight into the world of evaluation of stress-life fatigue using such results further illustrates the complexities of real-life engineering analysis, even though the presented case was simplified for the sake of clearness. Here is a list of the tips that arose throughout the text:
 
  * use and exercise your imagination
  * practise math
@@ -1105,29 +1129,28 @@ Back in college, we all learned how to solve engineering problems. And once we g
  * remember there are other methods beside finite elements
  * take into account that even within the finite element method, there is a wide variety of complexity in the problems that can be solved
 dnl * follow the “five whys rule” before compute anything, probably you do not need to
- * use engineering judgment and make sure understand Asimov’s [“wronger than wrong”](https://en.wikipedia.org/wiki/Wronger_than_wrong) concept
+dnl * use engineering judgement and make sure understand Asimov’s [“wronger than wrong”](https://en.wikipedia.org/wiki/Wronger_than_wrong) concept
  * play with your favourite FEM solver (mine is [CAEplex](https://caeplex.com)) solving simple cases, trying to predict the results and picturing the stress tensor and its eigenvectors in your imagination
  * prove (with pencil and paper) that even though the principal stresses are not linear with respect to summation, they are linear with respect to multiplication 
  * grab any stress distribution from any of your FEM projects, compute the iso-stress curves and the draw normal lines to them to get acquainted with SCLs
  * search online for “stress linearisation” (or “linearization” if you want) and then get a copy of ASME\ VIII Div\ 2 Annex 5-A
  * take into account that FEM solutions lead only to nodal equilibrium but not point-wise equilibrium
- * measure the time needed to generate grids of different sizes and kinds with your favourite mesher
- * learn this by heart: the complexity of a FEM problem is given mainly by the number of _nodes_, not by the number of elements
+dnl * measure the time needed to generate grids of different sizes and kinds with your favourite mesher
+dnl * learn this by heart: the complexity of a FEM problem is given mainly by the number of _nodes_, not by the number of elements
  * remember that welded materials with different thermal expansion coefficients may lead to fatigue under cyclic temperature changes
 dnl * if you have time, try to get out of your comfort zone and do more than what others expect from you (like parametric computations)
- * clone the [parametric tee repository](https://github.com/seamplex/tee), understand how the figures from\ [@sec:parametric] were built and expand them to cover “we might go on...” bullets
  * try to find an explanation of the results obtained, just like we did when we explained the parametric curves from\ [@fig:tee-MB ] with two opposing effects which were equal in magnitude around $d_b=5$\ inches
- * think thermal-mechanical plus earthquakes as “bake, break and shake” problems
+ * think thermal-mechanical plus earthquakes cases as “bake, break and shake” problems
  * understand why the elastic problem of the case study is still linear after all
  * keep in mind that finite-elements are a means to get an engineering solution, not an end by themselves
- * learn to write scripts to post-process FEM results (from a script-friendly FEM program)
+ * learn to write scripts to post-process FEM results (from an script-friendly open-source FEM program)
  * work under a [distributed version control](https://en.wikipedia.org/wiki/Distributed_version_control) system such as [Git](https://en.wikipedia.org/wiki/Git), even when just editing input files or writing reports
-dnl * clone the [environmental fatigue sample problem repository](https://bitbucket.org/seamplex/cufen) and obtain a nicely-formatted table with the results of the “EAF Sample Problem 2-Rev.\ 2 (10/21/2011)” from\ [@sec:in-air;@sec:in-water].
  * do not write ambiguous reports by replacing appropriate mathematical formulae with words just not to offend the illiterate
  * try to avoid [proprietary](https://en.wikipedia.org/wiki/Proprietary_software) programs and favour [free and open source](https://en.wikipedia.org/wiki/Free_and_open-source_software) ones.
 
 dnl You can ask for help in our mailing list at <wasora@seamplex.com>. There is a community of engineers willing to help you in case you get in trouble with the repositories, the script or the input files.
 
+
 About your favourite FEM program, ask yourself these two questions:
 
  1. Does your favourite FEM program’s manual say what the program does?
@@ -1137,5 +1160,32 @@ About your favourite FEM program, ask yourself these two questions:
 
 And finally, make sure that at the end of the journey from college theory to an actual engineering problem your conscience is clear knowing that there exists a report with your signature on it. That is why we all went to college in the first place.
 
+### Online stuff {#sec:online}
+
+Here is a list of sub-problems and stuff to play with.
+
+ * The pendulum-swing video from [@sec:intro]
+    - <https://youtu.be/Q-lKK4A2OzA>
+ * “On convergence of linearized stresses in an infinite pipe computed using the finite element method” from [@sec:infinite-pipe;@sec:linearity;@sec:infinite-pipe;@sec:infinite-pipe-fem]
+    - <https://www.seamplex.com/fino/doc/pipe-linearized/> 
+ * The three cubes from [@sec:linearity]
+    - case A: pure normal loads (<https://caeplex.com/p/d8fe>)
+    - case B: pure shear loads (<https://caeplex.com/p/b494>)
+    - case C: the combination of A & B (<https://caeplex.com/p/9899>)
+ * The animations of natural oscillations in [@fig:modes]
+    - <https://www.seamplex.com/docs/nafems4/mode1.webm>
+    - <https://www.seamplex.com/docs/nafems4/mode2.webm>
+    - ...
+ * The animations of natural oscillations in [@fig:case-mode]
+    - <https://www.seamplex.com/docs/nafems4/case-mode1.webm>
+    - <https://www.seamplex.com/docs/nafems4/case-mode2.webm>
+    - ...
+ * The parametric tee repository from [@sec:parametric]
+   - <https://github.com/seamplex/tee>
+ * The environmental fatigue sample problem repository from [@sec:in-air;@sec:in-water]
+   - <https://github.com/seamplex/cufen>
+
+See <https://www.seamplex.com/nafems> for new material, updated links and the full version of this case with many more details about the case and the associated mathematics.
+
 
 # Concluding Remarks

nafems4.md

@@ -11,16 +11,17 @@
 
 ## From college theory to an actual engineering problem {#sec:intro}
 
-First of all, please take this text as a written chat between you an me, i.e. an average engineer that has already taken the journey from college to performing actual engineering works using finite element analysis and has something to say about it. Picture yourself in a coffee bar, talking and discussing concepts and ideas with me. Maybe needing to go to a blackboard (or notepad?). Even using a tablet to illustrate some three-dimensional results. But always as a chat between colleagues.
+First of all, please take this text as a written chat between you an me, i.e. a practising engineer that has already taken the journey from college to performing actual engineering works using finite element analysis and has something to say about it. Picture yourself in a coffee bar, talking and discussing concepts and ideas with me. Maybe needing to go to a blackboard (or notepad?). Even using a tablet to illustrate some three-dimensional results. But always as a chat between colleagues.
 
-Please also note that I am not a mechanical engineer, although I shared many undergraduate courses with some of them. I am a nuclear engineer with a strong background in mathematics and computer programming. I went to college between 2002 and 2008. Probably a lot of things have changed since then---at least that is what these “millennial” guys and girls seem to be boasting about---but chances are we all studied solid mechanics and heat transfer with a teacher using a piece of chalk on a blackboard while we as students wrote down notes with pencils on paper sheets. And there is really not much that one can do with pencil and paper regarding mechanical analysis. Any actual case worth the time of an engineer needs to be more complex than an ideal canonical case with analytical solution.
+Please also note that I am not a mechanical engineer, although I shared many undergraduate courses with some of them. I am a nuclear engineer with a strong background in mathematics and computer programming. I went to college between 2002 and 2008. Probably a lot of things have changed since then---at least that is what these “millennial” guys and girls seem to be boasting about---but chances are we all studied solid mechanics and heat transfer with a teacher using a piece of chalk on a blackboard while we as students wrote down notes with pencils on paper sheets. And there is really not much that one can do with pencil and paper regarding mechanical analysis. Any actual case worth the time of an engineer needs to be more complex than an ideal canonical case with a closed-form solution.
 
-Whether you are a student or a seasoned engineer with many years of experience, you might recall from first year physics courses the introduction of the [simple pendulum](https://en.wikipedia.org/wiki/Pendulum) as a case study\ ([@fig:simple]). You learned that the period does not depend on the hanging mass because the weight and the inertia exactly cancelled each other. Also, that [Galileo](https://en.wikipedia.org/wiki/Galileo_Galilei) said (and [Newton proved](https://seamplex.com/wasora/doc/realbook/012-mechanics/)) that for small oscillations the period does not even depend on the amplitude. Someone showed you why it worked this way: because if\ $\sin \theta \approx \theta$ then the motion equations converge to an [harmonic oscillator](https://en.wikipedia.org/wiki/Harmonic_oscillator). It might have been a difficult subject for you back in those days when you were learning physics and calculus at the same time.
+Whether you are a student or a seasoned engineer with many years of experience, you might recall from first year physics courses the introduction of the [simple pendulum](https://en.wikipedia.org/wiki/Pendulum) as a case study\ ([@fig:simple]). You learned that the period does not depend on the hanging mass because the weight and the inertia exactly cancelled each other. Also, that [Galileo](https://en.wikipedia.org/wiki/Galileo_Galilei) said (and [Newton proved](https://seamplex.com/wasora/doc/realbook/012-mechanics/)) that for small oscillations the period does not even depend on the amplitude. Someone showed you why it worked this way: because if\ $\sin \theta \approx \theta$ (in radians) then the motion equations are identical to an [harmonic oscillator](https://en.wikipedia.org/wiki/Harmonic_oscillator). It might have been a difficult subject for you back in those days when you were learning physics and calculus at the same time.
 divert(-1)
 You might have later studied the [Lagrangian](https://en.wikipedia.org/wiki/Lagrangian_mechanics) and even the [Hamiltonian](https://en.wikipedia.org/wiki/Hamiltonian_mechanics) formulations, added a [parametric excitation](https://en.wikipedia.org/wiki/Parametric_oscillator) and analysed the [chaotic double pendulum](https://seamplex.com/wasora/doc/realbook/017-double-pendulum). divert(0)
-But it was probably after college, say when you took your first son to a swing on a windy day ([@fig:hamaca]), that you were faced with a real pendulum worth your full attention.
+But it was probably after college, say when you took your first child to a swing on a windy day ([@fig:hamaca]), that you were faced with a real pendulum worth your full attention.
 dnl Yes, this is my personal story but it could have easily been yours as well.
-The very same distance between what I imagined studying a pendulum was and what I saw that day at the swing (namely that the period _does_ depend on the hanging mass) is the same distance between the mechanical problems studied in college and the actual cases encountered during a professional engineer’s lifetime ([@fig:pipes]). In this regard, I am referring only to technical issues. The part of dealing with clients, colleagues, bosses, etc. which is definitely not taught in engineering schools (you can get a heads up in business schools, but again it would be a theoretical pendulum) is way beyond the scope of both this article and my own capacities.
+The very same difference between what I imagined studying a pendulum was and what I saw that day at the swing (namely that the period _does_ depend on the hanging mass) is the same difference between the mechanical problems studied in college and the actual cases encountered during a professional engineer’s lifetime ([@fig:pipes]). Real-life engineering problems are far more complex and involve far more information than what professors could write on a blackboard.
+In this regard, I am referring only to technical issues. The part of dealing with clients, colleagues, bosses, etc. which is definitely not taught in engineering schools (you can get an understanding in business schools, but again it would be a theoretical pendulum) is way beyond the scope of both this article and my own capacities.
 
 ::::: {#fig:pendulum}
 ![Simple pendulum](simple.svg){#fig:simple width=35%}\  
@@ -33,7 +34,7 @@ A [simple pendulum from college physics courses](https://www.seamplex.com/wasora
 ![College pipe](infinite-pipe.svg){#fig:infinite-pipe width=40%}
 ![Real-life pipe](isometric.svg){#fig:isometric width=58%}
 
-An infinitely-long pressurised thick pipe as taught in college and an isometric drawing of a section of a real-life piping system. The details of the isometric are not important individually but they are included to emphasize the fact that a real engineering drawing is far more complex and has far more information than what a professor can draw in a class blackboard.
+Left: what we are taught in college. Right: a real-life isometric drawing.
 :::::
 
 divert(-1)
@@ -44,21 +45,21 @@ Like the pendulums above, we will be swinging back and forth between a case stud
 Finite elements are like magic to me. I mean, I can follow the whole derivation of the equations, from the strong, weak and variational formulations of the equilibrium equations for the mechanical problem (or the energy conservation for heat transfer) down to the [algebraic multigrid](https://en.wikipedia.org/wiki/Multigrid_method) [preconditioner](https://en.wikipedia.org/wiki/Preconditioner) for the inversion of the [stiffness matrix](https://en.wikipedia.org/wiki/Stiffness_matrix) passing through [Sobolev spaces](https://en.wikipedia.org/wiki/Sobolev_space) and the [grid generation](https://en.wikipedia.org/wiki/Mesh_generation). Then I can sit down and program all these steps into a computer, including the [shape functions](https://en.wikipedia.org/wiki/Finite_element_method_in_structural_mechanics#Interpolation_or_shape_functions) and its derivatives, the assembly of the discretised stiffness matrix ([@sec:building]), the numerical solution of the system of equations ([@sec:solving]) and the computation of the gradient of the solution (@sec:stress-computation). Yet, the fact that all these a-priori unconnected steps give rise to pretty pictures that resemble reality is still astonishing to me.
 divert(0)
 
-Again, take all this information as coming from a fellow that has already taken such a journey from college’s pencil and paper to real engineering cases involving complex numerical calculations. And developing, in the meantime, both an actual working finite-element [back-end](https://www.seamplex.com/fino) and [front-end](https://www.caeplex.com) from scratch.dnl^[The dean of the engineering school I attended used to say “It is not the same to read than to write manuals, and we should aim at writing them.”]
+Again, take all this information as coming from a fellow that has already taken such a journey from college’s pencil and paper to real engineering cases involving complex numerical calculations. And developing, in the meantime, both an actual working finite-element [solver](https://www.seamplex.com/fino) and a [web-based pre and post-processor](https://www.caeplex.com) from scratch.dnl^[The dean of the engineering school I attended used to say “It is not the same to read than to write manuals, and we should aim at writing them.”]
 
 
-### Tips and tricks
+### Tips and hints
 
-There are some useful tricks that come handy when trying to solve a mechanical problem. Throughout this text, I will try to tell you some of them.
+There are some useful hints that come in handy when trying to solve a mechanical problem. Throughout this text, I will try to tell you some of them.
 
 
-One of the most important ones to use your _imagination_. You will need a lot of imagination to “see” what it is actually going on when analysing an engineering problem. This skill comes from my background in nuclear engineering where I had no choice but to imagine a [positron-electron annihilation](https://en.wikipedia.org/wiki/Electron%E2%80%93positron_annihilation) or an [spontaneous fission](https://en.wikipedia.org/wiki/Spontaneous_fission). But in mechanical engineering, it is likewise important to be able to imagine how the loads “press” one element with the other, how the material reacts depending on its properties, how the nodal displacements generate stresses (both normal and shear), how results converge, etc.
+One of the most important ones is to use your _imagination_. You will need a lot of imagination to “see” what it is actually going on when analysing an engineering problem. This skill comes from my background in nuclear engineering where I had no choice but to imagine a [positron-electron annihilation](https://en.wikipedia.org/wiki/Electron%E2%80%93positron_annihilation) or an [spontaneous fission](https://en.wikipedia.org/wiki/Spontaneous_fission). But in mechanical engineering, it is likewise important to be able to imagine how the loads “press” one element with the other, how the material reacts depending on its properties, how the nodal displacements generate stresses (both normal and shear), how results converge, etc.
 dnl And what these results actually mean besides the pretty-coloured figures.^[A former boss once told me “I need the CFD” when I handed in some results. I replied that I did not do computational fluid-dynamics but computed the neutron flux kinetics within a nuclear reactor core. He joked “I know, what I need are the _Colors For Directors_, those pretty-coloured figures along with your actual results to convince the managers.”]
 
 This journey will definitely need your imagination. We will peek a little bit into equations, numbers, plots, schematics, CAD geometries, 3D\ views, etc. Still, when the theory says “thermal expansion produces normal stresses” you have to picture in your head three little arrows pulling away from the same point in three directions, or whatever mental picture you have about what you understand thermally-induced stresses are. What comes to your mind when someone says that out of the nine elements of the [stress tensor](https://en.wikipedia.org/wiki/Cauchy_stress_tensor) ([@sec:tensor]) there are only six that are independent? Whatever it is, try to practice that kind of graphical thoughts with every new concept. Nevertheless, there will be particular locations of the text where imagination will be most useful. I will bring the subject up then and again.
 
 
-Another heads up is that we will be digging into some mathematics. Probably they would be simple and you would deal with them very easily. But chances are you do not like equations. No problem! Just ignore them for now. Read the text skipping them, it should work as well.
+Another point to observe is that we will be digging into some mathematics. Probably they would be simple and you would deal with them very easily. But chances are you do not like equations. No problem! Just ignore them for now. Read the text skipping them, it should work as well.
 Anyhow, here comes another experience tip: keep exercising mathematics. You have used [differences of squares](https://en.wikipedia.org/wiki/Difference_of_two_squares) in high school, didn’t you?. You know (or at least knew) how to [integrate by parts](https://en.wikipedia.org/wiki/Integration_by_parts). Do you remember what [Laplace transforms](https://en.wikipedia.org/wiki/Laplace_transform) are used for? Once in a while, perform a division of polynomials using [Ruffini’s rule](https://en.wikipedia.org/wiki/Ruffini's_rule). Or compute the second [derivative of the quotient of two functions](https://en.wikipedia.org/wiki/Quotient_rule). Whatever. It should be like doing crosswords on the newspaper. Grab those old physics college books and solve the exercises at the end of each chapter. All the effort  will pay off later on.
 
 One final comment: throughout the text I will be referring to “your favourite FEM program.” I bet you do have one. Mine is [CAEplex](https://caeplex.com) (it works on top of [Fino](https://www.seamplex.com/fino), which is free and open source). We will be using it to perform some tests and play a little bit. And we will also use it to think about what it means to use a FEM program to generate results that will eventually end up in a written project with your signature. Keep that in mind.
@@ -67,60 +68,60 @@ One final comment: throughout the text I will be referring to “your favourite
 ## Case study: reactors, pipes and fatigue {#sec:case}
 
 Piping systems in sensitive industries like nuclear or oil & gas should be designed and analysed following the recommendations of an appropriate set of codes and norms, such as the [ASME\ Boiler and Pressure Vessel Code](https://en.wikipedia.org/wiki/ASME_Boiler_and_Pressure_Vessel_Code).
-This code of practice was born during the late\ 19th century, before finite-element methods for solving partial differential equations were even developed. And much longer before they were available for the general engineering community. Therefore, much of the code assumes design and verification is not necessarily performed numerically but with paper and pencil (yes, like in college). However, it still provides genuine guidance in order to ensure pressurised systems behave safely and properly without needing to resort to computational tools. Combining finite-element analysis with the ASME code gives the cognisant engineer a unique combination of tools to tackle the problem of designing and/or verifying pressurised piping systems.
+This code of practice was born during the late\ 19th century, before finite-element methods for solving partial differential equations were even developed. And much longer before they were available for the general engineering community. Therefore, much of the code assumes design and verification is not necessarily performed numerically but with paper and pencil (yes, like in college). It provides guidance in order to ensure pressurised systems behave safely and properly without needing to resort to computational tools. Combining finite-element analysis with the ASME code gives the cognisant engineer a unique combination of tools to tackle the problem of designing and/or verifying pressurised piping systems.
 
-In the years following [Enrico Fermi](https://en.wikipedia.org/wiki/Enrico_Fermi)’s demonstration that a self-sustainable [fission reaction](https://en.wikipedia.org/wiki/Nuclear_fission) chain was possible (actually, in fact, after WWII was over), people started to build plants in order to transform the energy stored within the atoms nuclei into usable electrical power. They quickly reached the conclusion that high-pressure heat exchangers and turbines were needed. So they started to follow codes of practise like the aforementioned ASME\ B&PVC. They also realised that some requirements did not fit the needs of the nuclear industry. But instead of writing a new code from scratch, they added a new chapter to the existing body of knowledge: the celebrated [ASME Section\ III](https://en.wikipedia.org/wiki/ASME_Boiler_and_Pressure_Vessel_Code#ASME_BPVC_Section_III_-_Rules_for_Construction_of_Nuclear_Facility_Components).
+In the years following [Enrico Fermi](https://en.wikipedia.org/wiki/Enrico_Fermi)’s demonstration that a self-sustainable [fission reaction](https://en.wikipedia.org/wiki/Nuclear_fission) chain was possible (actually, in fact, after WWII was over), people started to build plants in order to transform the energy stored within an atom’s nucleus into usable electrical power. They quickly reached the conclusion that high-pressure heat exchangers and turbines were needed. So they started to follow codes of practise like the aforementioned ASME\ B&PVC. They also realised that some requirements did not fit the needs of the nuclear industry. But instead of writing a new code from scratch, they added a new chapter to the existing body of knowledge: the celebrated [ASME Section\ III](https://en.wikipedia.org/wiki/ASME_Boiler_and_Pressure_Vessel_Code#ASME_BPVC_Section_III_-_Rules_for_Construction_of_Nuclear_Facility_Components).
 
-After further years passed by, engineers (probably the same people that forked section\ III) noticed that fatigue in nuclear power plants was not exactly the same as in other piping systems. There were some environmental factors directly associated to the power plant that were not taken into account by the regular ASME code. Again, instead of writing a new code from scratch, people decided to add correction factors to the previously-amended body of knowledge. This is how (sometimes) knowledge evolves, and it is this kind of complexities that engineers are faced with during their professional lives. We have to admit it, it would be a very difficult task to re-write everything from scratch every time something changes.
+After further years passed by, engineers (probably the same people that wrote section\ III) noticed that fatigue in nuclear power plants was not exactly the same as in other piping systems. There were some environmental factors directly associated to the power plant that were not taken into account by the regular ASME code. Again, instead of writing a new code from scratch, people decided to add correction factors to the previously-amended body of knowledge. This is how (sometimes) knowledge evolves, and it is these kinds of complexities that engineers are faced with during their professional lives. We have to admit it, it would be a very difficult task to re-write everything from scratch every time something changes.
 
-![Three-dimensional CAD model for the piping system of\ [@fig:isometric]](real-piping.png){#fig:real-life width=75%}
+![Three-dimensional CAD model for the piping system of\ [@fig:isometric]](real-piping.png){#fig:real-life width=100%}
 
 
 ### Nuclear reactors
 
-In each of the countries that have at least one nuclear power plant there exists a national regulatory body who is responsible for allowing the owner to operate the reactor. These operating licenses are time-limited, with a range that can vary from 25 to 60 years, depending on the design and technology of the reactor. Once expired, the owner might be entitled to an extension, which the regulatory authority can accept provided it can be shown that a certain (and very detailed) set of safety criteria are met. One particular example of requirements is that of fatigue in pipes, especially those that belong to systems that are directly related to the reactor safety.
+In each of the countries that have at least one nuclear power plant there exists a national regulatory body who is responsible for allowing the owner to operate the reactor. These operating licenses are time-limited, with a range that can vary from 25 to 60 years, depending on the design and technology of the reactor. Once expired, the owner might be entitled to an extension, which the regulatory authority can grant provided it can be shown that a certain (and very detailed) set of safety criteria are met. One particular example of requirements is that for fatigue in pipes, especially those that belong to systems that are directly related to the reactor safety.
 
 ### Pressurised pipes
 
-How come that pipes are subject to fatigue? Well, on the one hand and without getting into many technical details, the most common nuclear reactor design uses liquid water as coolant and moderator. On the other hand, nuclear power plants cannot by-pass the thermodynamics of the [Carnot cycle](https://en.wikipedia.org/wiki/Carnot_cycle), and in order to maximise the efficiency of the conversion of the energy stored in the uranium nuclei into electricity we need to reach temperatures as high as possible. So, if we want to have liquid water in the core as hot as possible, we need to increase the pressure. The limiting temperature and pressure are given by the [critical point of water](https://en.wikipedia.org/wiki/Critical_point_(thermodynamics)), which is around 374ºC and 22\ MPa. It is therefore expected to have temperatures and pressures near those values in many systems of the plant, especially in the primary circuit and those that directly interact with it, such as pressure and inventory control system, decay power removal system, feedwater supply system, emergency core-cooling system, etc.
+Why are pipes are subject to fatigue? Well, on the one hand and without getting into many technical details, the most common nuclear reactor design uses liquid water as coolant and moderator. On the other hand, nuclear power plants cannot by-pass the thermodynamics of the [Carnot cycle](https://en.wikipedia.org/wiki/Carnot_cycle), and in order to maximise the efficiency of the conversion of the energy stored in the uranium nuclei into electricity we need to reach temperatures as high as possible. So, if we want to have liquid water in the core as hot as possible, we need to increase the pressure. The limiting temperature and pressure are given by the [critical point of water](https://en.wikipedia.org/wiki/Critical_point_(thermodynamics)), which is around 374ºC and 22\ MPa. It is therefore expected to have temperatures and pressures near those values in many systems of the plant, especially in the primary circuit and those that directly interact with it, such as pressure and inventory control system, decay power removal system, feedwater supply system, emergency core-cooling system, etc.
 
-Nuclear power plants are not always working at 100% of their maximum power capacity. They need to be maintained and refuelled, they may undergo operational (and some incidental) transients, they might operate at a lower power due to load following conditions, etc. These transient cases involve changes both in temperatures and in pressures that the pipes are subject to, which in turn give rise to changes in the stresses within the pipes. As the transients are postulated to occur cyclically during a number of times during the life-time of the plant (plus its extension period), mechanical fatigue in these piping systems may arise, especially at the interfaces between materials with different thermal expansion coefficients.
+Nuclear power plants are not always working at 100% of their maximum power capacity. They need to be maintained and refuelled, they may undergo operational (and some incidental) transients, they might operate at a lower power due to load following conditions, etc. These transient cases involve temperatures and pressures changes in the pipes, which in turn give rise to changes in the stresses within the pipes. As the transients are postulated to occur cyclically during a number of times throughout the life-time of the plant (plus its extension period), mechanical fatigue in these piping systems may arise, especially at the interfaces between dissimilar materials e.g. with different thermal expansion coefficients.
 
 An important part of the analysis that almost always applies to nuclear power plants but usually also to other installations is the consideration of a possible seismic event. Given a postulated design earthquake, both the civil structures and the piping system itself need to be able to withstand such a load, even if it occurs at the moment of highest mechanical demand during one of the operational transients.
 
 ### Fatigue {#sec:fatigue}
 
-Mechanical systems can fail due to a wide variety of reasons. The effect known as fatigue can create,  migrate and grow microscopic cracks at the atomic level, called [dislocations](https://en.wikipedia.org/wiki/Dislocations). Once these cracks reach a critical size, then the material fails catastrophically even under stresses much lower than the tensile strength limits. There are not complete mechanistic models from first principles which can be used in general situations, and those that exist are very complex and hard to use. There are two main ways to approach practical fatigue assessment problems using experimental data very much like the [Hooke’s Law](https://en.wikipedia.org/wiki/Hooke's_law) experiment:
+Mechanical systems can fail for a wide variety of reasons. The effect known as fatigue can initiate and grow microscopic cracks at the grain level, called [dislocations](https://en.wikipedia.org/wiki/Dislocations). Cracks initiation grow at levels well below yield stresses. Once these cracks reach a critical size, then the material fails catastrophically. There are currently not complete mechanical models describing fatigue from first principles. Therefore an empirical approach is used. There are two main ways to approach practical fatigue assessment problems using experimental data:
 
  1. stress life, or
  2. strain life
 
 The first one is suitable for cases where the loads are nowhere near the yield stress of the material. When plastic deformation is expected to occur, strain-life methods ought to be employed.
 
-For the case study, as the loads come principally from operational loads, the ASME\ stress-life approach should be used. The stress amplitude of a periodic cycle can be related to the number of cycles where failure by fatigue is expected to occur. For each material, this dependence can be computed using normalised tests and a family of “fatigue curves” like the one depicted in\ [@fig:SN] (also called $S$-$N$ curve) for different temperatures can be obtained.
+For the case study, as the loads come principally from operational loads, the ASME\ stress-life approach should be used. The stress amplitude (AS) of a periodic cycle can be related to the number of cycles ($N$) where failure by fatigue is expected to occur. For each material, this dependence can be computed using standardised tests and a family of “fatigue curves” like the one depicted in\ [@fig:SN] (also called $S$-$N$ curve) for different temperatures can be obtained.
 
 
 ![A fatigue or $S$-$N$ curve for two steels.](SN.svg){#fig:SN width=75%}
 
-It should be noted that the fatigue curves are obtained in a particular load case, namely purely-periodic and one-dimensional, which cannot be directly generalised to other three-dimensional cases. Also, any real-life case will be subject to a mixture of complex cycles given by a stress time history and not to pure periodic conditions. The application of the curve data implies a set of simplifications and assumptions that are translated into different possible “rules” for composing real-life cycles. There also exist two safety factors which increase the stress amplitude and reduce the number of cycles respectively. All these intermediate steps render the analysis of fatigue into a conservative computation scheme. Therefore, when a fatigue assessment performed using the fatigue curve method arrives at the conclusion that “fatigue is expected to occur after ten thousand cycles” what it actually means is “we are sure fatigue will not occur before ten thousand cycles, yet it may not occur before one hundred thousand or even more.”  
+It should be noted that the fatigue curves are obtained in a particular load case, namely purely-periodic and one-dimensional, which cannot be directly generalised to other three-dimensional cases. Also, any real-life case will be subject to a mixture of complex cycles given by a stress time history and not to pure periodic conditions. The application of the curve data implies a set of simplifications and assumptions that are translated into different possible “rules” for composing real-life cycles. The ASME $S$-$N$ curve also adopt tow safety factors which increase the stress amplitude and reduce the number of cycles respectively. All these intermediate steps render the analysis of fatigue into a conservative computation scheme. Therefore, when a fatigue assessment performed using the fatigue curve method arrives at the conclusion that “fatigue is expected to occur after ten thousand cycles” what it actually means is “we are sure fatigue will not occur before ten thousand cycles, yet it may not occur before one hundred thousand or even more.”  
 
 ## Solid mechanics, or what we are taught at college
 
-So, let us start our journey. Our starting place: undergraduate solid mechanics courses. Our goal: to obtain the internal state of a solid subject to a set of movement restrictions and loads (i.e. to solve the solid mechanics problem). Our first step: [Newton’s laws of motion](https://en.wikipedia.org/wiki/Newton's_laws_of_motion). For our purposes, we can recall them here like this:
+So, let us start our journey. Our starting place: undergraduate solid mechanics courses. Our goal: to obtain the internal state of a solid subject to a set of displacement restrictions and loads (i.e. to solve the solid mechanics problem). 
+It was [Augustin-Louis Cauchy](https://en.wikipedia.org/wiki/Augustin-Louis_Cauchy) who formulated for the first time the elasticity equations  we use today. We need to simultaneously solve
 
- 1. a solid is in equilibrium if it is not moving in at least one inertial coordinate system,
- 2. in order for a solid not to move, the sum of all the forces ought to be equal to zero, and
- 3. for every external load there exists an internal reaction with the same magnitude but opposite direction.
+ 1. equilibrium equations: the fact that external loads need to be balanced by internal stresses,
+ 2. constitutive equations: how strains generate stresses within the solid’s material(s), and
+ 3. compatibility equations: make sure that the deformed body satisfies the displacement restrictions.
 
-We have to accept that there is certain intellectual beauty when complex stuff can be expressed in such simple terms. Yet, from now on, everything can be complicated at will. We can take the mathematical path like [D’Alembert](https://en.wikipedia.org/wiki/Jean_le_Rond_d'Alembert) and his virtual displacements ideas (in his mechanical treatise, he brags that he does not need to use a single figure throughout the book). Or we can go graphical following [Culmann](https://en.wikipedia.org/wiki/Carl_Culmann). Or whatever other logic reasoning to end up with a set of actual equations which we need to solve in order to obtain engineering results.
+ 
+### Cauchy’s stress tensor {#sec:tensor}
 
-### The stress tensor {#sec:tensor}
-
-In any case, what we should understand (and imagine) is that external forces lead to internal stresses. And in any three-dimensional body subject to such external loads, the best way to represent internal stresses is through a $3 \times 3$ _stress tensor_. This is the first point in which we should not fear mathematics. Trust me, it will pay back later on.
+In any case, what we need to understand (and imagine) from point\ #1 above is that external forces lead to internal stresses. And in any three-dimensional body subject to such external loads, the best way to represent internal stresses is through a $3 \times 3$ _stress tensor_. This is the first point in which we should not fear mathematics. Trust me, let’s follow good old Augustin-Louis. It will pay back later on.
 
 Does the term _tensor_ scare you? It should not. A [tensor](https://en.wikipedia.org/wiki/Tensor) is a general mathematical object and might get complex when dealing with many dimensions (as those encountered in weird stuff like [string theory](https://en.wikipedia.org/wiki/String_theory)), but we will stick here to second-order tensors. They are slightly more complex than a vector, and I assume you are not afraid of vectors, are you? If you recall fresh-year algebra courses, a [vector](https://en.wikipedia.org/wiki/Vector_(mathematics_and_physics)) somehow generalises the idea of a [scalar](https://en.wikipedia.org/wiki/Scalar_(mathematics)) in the following sense: a given vector $\mathbf{v}$ can be projected into any direction $\mathbf{n}$ to obtain a scalar $p$. We call this scalar $p$ the “projection” of the vector $\mathbf{v}$ in the direction $\mathbf{n}$. Well, a tensor can be also projected into any direction $\mathbf{n}$. The difference is that instead of a scalar, a vector is now obtained.
 
-Let me introduce then the three-dimensional [stress tensor](https://en.wikipedia.org/wiki/Cauchy_stress_tensor) (a.k.a [Cauchy](https://en.wikipedia.org/wiki/Augustin-Louis_Cauchy) tensor):
+So, let me introduce then the three-dimensional [Cauchy’s stress tensor](https://en.wikipedia.org/wiki/Cauchy_stress_tensor):
 
 $$
 \begin{bmatrix}
@@ -130,13 +131,13 @@ $$
 \end{bmatrix}
 $$
 
-It looks (and works) like a regular $3 \times 3$ matrix. Some brief comments about it:
+It looks (and works) like a regular $3 \times 3$ matrix. Indeed, we will take advantage of this matrix-like behaviour in [@sec:linearity] below. Some brief comments about it:
 
- * The $\sigma$s are normal stresses, i.e. they try to stretch or tighten the material.
+ * The $\sigma$s are normal stresses, i.e. they try to stretch or tighten the material (depending on the sign).
  * The $\tau$s are shear stresses, i.e. they try to twist the material.
  * Due to rotational equilibrium requirements the conjugate shear stresses should be equal: $\tau_{xy} = \tau_{yx}$, $\tau_{yz} = \tau_{zy}$,  and $\tau_{zx} = \tau_{xz}$. Therefore, the stress tensor is symmetric i.e. there are only six independent elements.
  * The elements of the tensor depend on the orientation of the coordinate system.
- * There exists a particular coordinate system in which the stress tensor is diagonal, i.e. all the shear stresses are zero. In this case, the three diagonal elements are called the [principal stresses](https://en.wikipedia.org/wiki/Principal_stresses), which happen to be the three [eigenvalues](https://en.wikipedia.org/wiki/Eigenvalues_and_eigenvectors) of the stress tensor. The basis of the coordinate system that renders the tensor diagonal are the eigenvectors.
+ * There exists a particular coordinate system in which the stress tensor is diagonal, i.e. where all the shear stresses are zero. In this case, the three diagonal elements are called the [principal stresses](https://en.wikipedia.org/wiki/Principal_stresses), which happen to be the three [eigenvalues](https://en.wikipedia.org/wiki/Eigenvalues_and_eigenvectors) of the stress tensor. The basis of the coordinate system that renders the tensor diagonal are the eigenvectors.
 
 
 What does this all have to do with mechanical engineering? Well, once we know what the stress tensor is for every point of a solid, in order to obtain the internal forces per unit area acting in a plane passing through that point and with a normal given by the direction $\mathbf{n}$, all we have to do is “project” the stress tensor through $\mathbf{n}$. In plain simple words:
@@ -146,9 +147,46 @@ What does this all have to do with mechanical engineering? Well, once we know wh
 
 ### An infinitely-long pressurised pipe {#sec:infinite-pipe}
 
-Let us proceed to our second step, and consider the infinite pipe subject to uniform internal pressure already introduced in\ [@fig:infinite-pipe]. Actually, we are going to solve the mechanical problem on an infinite hollow cylinder, which looks like pipe. This case is usually tackled in college courses, and chances are you already solved it. In fact, the first (and simpler) problem is the “thin cylinder problem.” Then, the “thick cylinder problem” is introduced (the one we solve below), which is slightly more complex. Nevertheless, it has an analytical solution.^[A detailed analysis of the analytical solution and how results obtained with finite elements converge with respect to the mesh size can be found in the report “On convergence of linearized stresses in an infinite pipe computed using the finite element method.” See [@sec:online].] For the present case, let us consider an infinite pipe (i.e. a hollow cylinder) of internal radius $a$ and external radius $b$ with uniform mechanical properties---Young’s modulus $E$ and Poisson’s ratio $\nu$---subject to an internal uniform pressure $p$.
+Let us proceed to our second step, and consider the infinite pipe subject to uniform internal pressure already introduced in\ [@fig:infinite-pipe]. Actually, we are going to solve the mechanical problem on an infinite hollow cylinder, which looks like a pipe. This case is usually tackled in college courses, and chances are you already solved it. In fact, the first (and simpler) problem is the “thin cylinder problem.” Then, the “thick cylinder problem” is introduced (the one we solve below), which is slightly more complex. Nevertheless, it has an analytical solution. For the present case, let us consider an infinite pipe (i.e. a hollow cylinder) of internal radius $a$ and external radius $b$ with uniform mechanical properties---Young’s modulus $E$ and Poisson’s ratio $\nu$---subject to an internal uniform pressure $p>0$. [@Fig:timoshenko-cyl] provides a classical illustration of this problem. What follows is more or less what we are taught in school: some equations with a brief explanation of the results. And then we move on to the next subject.
+
+::::: {#fig:timoshenko}
+![Thick cylinder under pressure](fig41.svg){#fig:timoshenko-cyl width=48%}
+![Equilibrium in cylindrical coordinates.](fig40.svg){#fig:timoshenko-eq width=48%}
+
+Figures from [Semyon Timoshenko](https://en.wikipedia.org/wiki/Semyon_Timoshenko) classical book “Theory of elasticity.”
+:::::
+
+Given the cylindrical symmetry of the problem, there is only one independent variable, namely the radial coordinate $r$. Moreover, there are only two displacement fields that need to be considered: the axial $u_a(r)$ and the radial $u_r(r)$. The former is identically zero due to the fact that the cylinder is infinite ([plane strain condition](https://www.ramsay-maunder.co.uk/downloads/NBR03.pdf)). The equilibrium equation along the radial direction $r$, also known as the [Lamé](https://en.wikipedia.org/wiki/Gabriel_Lam%C3%A9) equation, can be derived with the aid of [@fig:timoshenko-eq] as
+
+$$ \frac{d\sigma_r}{dr} + \frac{\sigma_r(r) - \sigma_\theta(r)}{r} = 0 $${#eq:equilibrium}
+
+Defining the strains as $\epsilon_r= du_r/dr$ and $\epsilon_\theta = u/r$, the constitutive equations of an isotropic linear material are
 
-What follows is more or less what we are taught in school: some equations with a brief explanation of the results. And then we move on to the next subject.
+\begin{align*}
+\sigma_r &= \frac{E}{(1-\nu)(1-2\nu)} \cdot \Big[ (1-\nu) \cdot \epsilon_r + \nu \cdot \epsilon_\theta \Big] \\
+\sigma_\theta &= \frac{E}{(1-\nu)(1-2\nu)} \cdot  \Big[ \nu \cdot \epsilon_r + (1-\nu) \cdot \epsilon_r \Big] \\
+\end{align*}
+then the differential equation [-@eq:equilibrium] can be casted in terms of the axial displacement $u_r(r)$ as
+
+$$
+\frac{d^2 u}{dr^2} + \frac{1}{r} \cdot \frac{du}{dr} - \frac{u}{r^2} = 0
+$$
+that has the general solution
+
+$$
+u(r) = c_1 \cdot r + \frac{c_2}{r}
+$$
+
+For the boundary conditions of the particular problem that the radial stress should be equal to the negative of the internal pressure at $r=a$ and null at $r=b$, the axial displacement has the particular solution
+
+$$
+u_r(r) = p \cdot \frac{1+\nu}{E} \cdot \frac{a^2}{b^2-a^2} \cdot \left[ 1-2\nu + \frac{b^2}{r^2} \right]\cdot r
+$${#eq:ur}
+
+Exercises for the reader:
+
+ 1. Think (imagine) and answer why the radial stress should be equal to the negative of the internal pressure at $r=a$ and not the actual pressure\ $p$. Hint: it does not have anything to do with actions and reactions.
+ 2. Verify that the radial displacement $u_r(r)$ satisfies the equilibrium and constitutive equations.
 
 
 #### Displacements
@@ -157,13 +195,9 @@ Remember that when any solid body is subject to external forces, it has to react
 
  * There are no longitudinal displacements\ $u_l$ because the pipe is infinite in the axial direction.
  * There are no azimuthal displacements\ $u_\theta$ because the pipe is fully symmetric around the axis.
- * There are only radial displacements\ $u_r$ and they depend only on the radial coordinate\ $r$ and not on the axial position\ $z$ or on the azimuthal angle\ $\theta$. These displacements are
+ * There are only radial displacements\ $u_r$ and they depend only on the radial coordinate\ $r$ and not on the axial position\ $z$ or on the azimuthal angle\ $\theta$. 
 
- $$
- u_r(r) = p \cdot \frac{1+\nu}{E} \cdot \frac{a^2}{b^2-a^2} \cdot \left[ 1-2\nu + \frac{b^2}{r^2} \right]\cdot r
- $$
-
-What does this mean? Well, that overall the whole pipe expands a little bit radially with the inner face being displaced more than the external surface (use your imagination!). How much?
+What does this mean? Well, that overall the whole pipe expands a little bit radially with the inner face being displaced more than the external surface (use your imagination or wait until we get to [@fig:ur]). How does it expand?
 
  1. linearly with the pressure, i.e. twice the pressure means twice the displacement, and
  2. inversely proportional to the Young’s Modulus\ $E$ divided by $1+\nu$, i.e. the more resistant the material, the less radial displacements.
@@ -234,12 +268,18 @@ Discretisation of a spatial domain, For the same number of cells, unstructured g
 Back to the three numerical methods, we must say that finite differences is based on approximating derivatives (i.e. differentials) by incremental quotients (i.e. differences). The second one heavily relies on geometrical ideas rather than on pure mathematical grounds. Finally, our beloved finite elements are more strict in the mathematical sense. Actually, a complete derivation of the finite element method can be written in a textbook without requiring a single figure, just like D’Alembert did more than two centuries ago. In any case, it is important to note that finite differences and elements compute results at the _nodes_ of a mesh, whilst finite volumes compute results at the _cells_ of a mesh. Finally, any method may be used in structured grids but only finite elements and volumes are especially suited for working with unstructured grids.
 
 There are technical reasons that justify why the finite element method is the king of mechanical analysis. But that does not mean that other methods may be employed. For instance, fluid mechanics are generally better solved using finite volumes. And further other combinations may be found in the literature.
-
 divert(0)
 
+[Fig:ur] shows how some results obtained with the finite element method using different number and order of elements compare to the analytical solution from the previous section. A detailed analysis of the analytical solution and how results obtained with finite elements converge with respect to the mesh size can be found in the report “On convergence of linearised stresses in an infinite pipe computed using the finite element method.” See [@sec:online].
+
+![Comparison of $u_r(t)$ computed with finite elements and\ [@eq:ur].](ur.svg){#fig:ur width=90%}
+
+
 Before proceeding, I would like to make two comments about common nomenclature. The first one is that if we exchanged the words “volumes” and “elements” in all the written books and articles, nobody would notice the difference. There is nothing particular in both theories that can justify why FVM use “volumes” and FEM use “elements”. Actually volumes and elements are the same geometric constructions. As far as I know, the names were randomly assigned.
 
-The second one is more philosophical and refers to the word “simulation” which is often used to refer to solving a problem using a numerical scheme such as the finite element method. [I am against at using this word for this endeavour](https://www.seamplex.com/blog/say-modeling-not-simulation.html). The term simulation has a connotation of both “pretending” and “faking” something, that is definitely not what we are doing when we solve an engineering problem with finite elements. Sure, there are some cases in which we simulate, such as using the [Monte Carlo method](https://en.wikipedia.org/wiki/Monte_Carlo_method) (originally used by Fermi as an attempt to understand how neutrons behave in the core of nuclear reactors). But when solving deterministic mechanical engineering problems I would rather say “modelling” than “simulation.”
+
+
+dnl The second one is more philosophical and refers to the word “simulation” which is often used to refer to solving a problem using a numerical scheme such as the finite element method. [I am against at using this word for this endeavour](https://www.seamplex.com/blog/say-modeling-not-simulation.html). The term simulation has a connotation of both “pretending” and “faking” something, that is definitely not what we are doing when we solve an engineering problem with finite elements. Sure, there are some cases in which we simulate, such as using the [Monte Carlo method](https://en.wikipedia.org/wiki/Monte_Carlo_method) (originally used by Fermi as an attempt to understand how neutrons behave in the core of nuclear reactors). But when solving deterministic mechanical engineering problems I would rather say “modelling” than “simulation.”
 
 divert(-1)
 ### Kinds of finite elements {#sec:kinds}
@@ -349,7 +389,7 @@ CAD model of a piping system with a 3/4-inch weldolet-type fork (stainless steel
 :::::
 
 
-![Location of six SCLs defined to analyse fatigue around a junction.](weldolet-scls-n.png){#fig:weldolet-scls width=75%}
+![Location of six SCLs defined to analyse fatigue around a junction.](weldolet-scls-n.png){#fig:weldolet-scls width=90%}
 
 
 On the one hand, a reasonable number of nodes (it is the number of nodes that defines the problem size, not the number of elements) in order to get a decent grid is around 200k for each system. On the other hand, solving a couple of dozens of transient heat transfer problems (which we cannot avoid due to the large thermal inertia of the pipes) during a few thousands of seconds over a couple hundred of thousands of nodes might take more time and storage space to hold the results than we might expect.
@@ -362,11 +402,11 @@ We can then merge this idea by Asimov with an adapted version of the [Saint-Vena
  2. assume the part of the full system which is not contained in the reduced mesh is at an uniform (though not constant) temperature equal to the average of the inner and outer temperatures at each side of the reduced mesh.
 
 
-![An example case where the SCLs are located around the junction between stainless-steel valves and carbon steel pipes at both sides of the material interface in the vertical plane both at the top and at the bottom of the pipe.](valve-scls1-n.png){#fig:valve width=65%}
+![An example case where the SCLs are located around the junction between stainless-steel valves and carbon steel pipes at both sides of the material interface in the vertical plane both at the top and at the bottom of the pipe.](valve-scls1-n.png){#fig:valve width=85%}
 
 As an example, let us consider the system depicted in\ [@fig:valve] where there is a stainless-carbon steel interface at the discharge of the valves. Instead of solving the transient heat-conduction problem with the internal temperature of the pipes equal to the temperature of the water in the reference transient condition of the power plant and an external condition of natural convection to the ambient temperature in the whole mesh, a reduced model consisting of half of one of the two valves and a small length of the pipes at both the valve inlet and outlet is used. Once the temperature distribution\ $\hat{T}(\mathbf{x},t)$ for each time is obtained in the reduced mesh ([@fig:valve-temp], which has the origin at the centre of the valve), the actual temperature distribution\ $T(\mathbf{x},t)$ is computed by an algebraic generalisation of $\hat{T}(\mathbf{x},t)$ in the full coordinate system. As stated above, those locations which are not covered by the reduced model are generalised with a time-dependent uniform temperature which is the average of the inner and outer temperatures at the inlet and outlet of the reduced mesh.
 
-![Mesh of half-a-symmetric-valve refined around the interface where the transient heat conduction problem is solved.](valve-temp-commented.png){#fig:valve-temp width=80%}
+![Mesh of half-a-symmetric-valve refined around the interface where the transient heat conduction problem is solved.](valve-temp-commented.png){#fig:valve-temp width=100%}
 
 
 Note that there is no need to have a one-to-one correspondence between the elements from the reduced mesh with the elements from the original one. Actually, the reduced mesh contains first-order elements whilst the former has second-order elements. Also the grid density is different. Nevertheless, the finite-element solver [Fino](https://www.seamplex.com/fino) used to solve both the heat and the mechanical problems, allows to read functions of space and time defined over one mesh and continuously evaluate and use them into another one even if the two grids have different elements, orders or even dimensions. In effect, in the system from\ [@fig:real-life] the material interface is between a orifice plate made in stainless steel that is welded to a carbon-steel pipe ([@fig:real]). The thermal problem can be modelled using a two-dimensional axi-symmetric grid and then generalised to the full three-dimensional mesh using the algebraic manipulation capabilities provided by [Fino](https://www.seamplex.com/fino) (actually by [wasora](https://www.seamplex.com/wasora)) as shown in\ [@fig:real].
@@ -374,7 +414,7 @@ Note that there is no need to have a one-to-one correspondence between the eleme
 dnl ![Full 3D mesh, carbon steel is magenta, stainless steel is green.](real-mesh2.png){#fig:real-mesh2 width=70%}
 dnl ![Temperature distribution for a certain time within the transient computed on a reduced two-dimensional axi-symmetric mesh modelling half the orifice plate and a length of the carbon pipe.](real-temp.png){#fig:real-temp width=34%}
 
-![Temperature distribution for a certain time within the transient computed on a reduced two-dimensional axi-symmetric mesh modelling half the orifice plate and a length of the carbon pipe and generalisation of the temperature to the full three-dimensional grid.](real-gen.png){#fig:real width=85%}
+![Temperature distribution for a certain time within the transient computed on a reduced two-dimensional axi-symmetric mesh modelling half the orifice plate and a length of the carbon pipe and generalisation of the temperature to the full three-dimensional grid.](real-gen.png){#fig:real width=100%}
 
 
 ### Seismic loads {#sec:seismic}
@@ -545,7 +585,7 @@ ans =
    5.767255
 ```
 
-Did I convince you? More or less, right? The third eigenvalue seems to fit. Let us not throw all of our beloved linearity away and dig in further into the subject. There are still two important issues to discuss which can be easily addressed using fresh-year linear algebra (remember not to fear maths!). First of all, even though principal stresses are not linear with respect to the sum they are linear with respect to pure multiplication. Once more, think what happens to the the eigenvalues and eigenvectors of a single stress tensor as all its elements are scaled up or down by a real scalar. They are the same! So, for example, the [Von\ Mises stress](https://en.wikipedia.org/wiki/Von_Mises_yield_criterion) (which is a combination of the principal stresses) of a beam loaded with a force\ $\alpha \cdot \mathbf{F}$ is\ $\alpha$ times the stress of the beam loaded with a force\ $\mathbf{F}$. Please test this hypothesis by playing with your favourite FEM solver. Or even better, take a look at the stress invariants $I_1$, $I_2$ and $I_3$ (you can search online or peek into the source code of [Fino](https://www.seamplex.com/fino), grep for the routine called `fino_compute_principal_stress()`) and see (using paper and pencil!) how they scale up if the individual elements of the stress tensor are scaled by a real factor\ $\alpha$.
+<!-- Did I convince you? More or less, right? The third eigenvalue seems to fit. Let us not throw all of our beloved linearity away and dig in further into the subject. There are still two important issues to discuss which can be easily addressed using fresh-year linear algebra (remember not to fear maths!). First of all, even though principal stresses are not linear with respect to the sum they are linear with respect to pure multiplication. Once more, think what happens to the the eigenvalues and eigenvectors of a single stress tensor as all its elements are scaled up or down by a real scalar. They are the same! So, for example, the [Von\ Mises stress](https://en.wikipedia.org/wiki/Von_Mises_yield_criterion) (which is a combination of the principal stresses) of a beam loaded with a force\ $\alpha \cdot \mathbf{F}$ is\ $\alpha$ times the stress of the beam loaded with a force\ $\mathbf{F}$. Please test this hypothesis by playing with your favourite FEM solver. Or even better, take a look at the stress invariants $I_1$, $I_2$ and $I_3$ (you can search online or peek into the source code of [Fino](https://www.seamplex.com/fino), grep for the routine called `fino_compute_principal_stress()`) and see (using paper and pencil!) how they scale up if the individual elements of the stress tensor are scaled by a real factor\ $\alpha$. -->
 
 The other issue is that even though in general the eigenvalues of the sum of two matrices are not the same as the eigenvalues of the matrix sum, there are some cases when they are. In effect, if two matrices\ $A$ and\ $B$ commute, i.e. their product is commutative
 
@@ -1045,9 +1085,6 @@ A design-basis earthquake was assumed to occur briefly during one second (sic) a
 
 To compute the global usage factor, we first need to find all the combinations of local extrema pairs and then sort them in decreasing order of stress difference. For example, the largest stress amplitude is found between $i=447$ and $i=694$. The second one is 447--699. Then 699--1020, 699--899, etc. For each of these pairs, defined by the indexes\ $i_{1,j}$ and $i_{2,j}$, a partial usage factor\ $U_j$ should computed. The stress amplitude\ $S_{\text{alt},j}$ which should be used to enter into the $S$-$N$ curve is
 
-\pagebreak
-
-
 $$
 S_{\text{alt},j} = \frac{1}{2} \cdot k_{e,j} \cdot \left| MB^\prime_{i_{1,j}} - MB^\prime_{i_{2,j}} \right| \cdot \frac{E_\text{SN}}{E(T_{\text{max}_j})}
 $$