*6.2. Privacy Preserving*

In the proposed scheme, if we consider the encrypted text by a home appliance

$$c\_{ip} = \left[1 + n \cdot \sigma\_{\dot{\jmath}} \cdot \left(\mathbf{x}\_{i} \cdot \sigma\_{0} + \mathbf{x}\_{i}^{2}\right)\right] \cdot H(T\_{s})^{n \cdot p\_{i}} \mod n^{2} \tag{30}$$

and that of the aggregated cipher text [82]

$$C\_p^\* = \left(1 + n \cdot \sum\_{\substack{i=1, i \neq a}}^N \sigma\_j^\* \cdot (x\_i \cdot \sigma\_0 + x\_i^2) \right) \cdot \prod\_{\substack{i=1, i \neq a}}^{N+1} H(T\_s)^{n \cdot p\_i} \bmod n^2$$

and take *σj* · (*xi* · *σ*0+*x*2 *i* ) as a *message* and ·*H*(*Ts*)·*pi* as a random number *rand* and ∑ *N i*=1,*i*=*a σ*∗ *j* · (*xi* · *σ*0 + *x*2 *i* ) as *Message* and ∏ *N*+1 *i*=1,*i*=*a H*(*Ts*)·*pi* as *Rand* then

$$x\_{ip} = \left[1 + n \cdot \overline{\text{message}\,\theta}\right] \cdot \overline{rand}^n \text{ mod } n^2 \tag{31}$$

and

$$\mathbb{C}\_p^\* = \left(1 + n \cdot \overbrace{\mathbb{M}essage}\right) \cdot \overline{\mathbb{R}and}^n \text{ mod } n^2$$

are valid Paillier ciphers. Under the chosen plaintext attack, Paillier encryption is indistinguishable, and an external attacker cannot obtain the exact message. The sink node may be curious about the exact message. However, without the knowledge of secret keys *p*0 and *λ*, the sink node has no knowledge of the information. Whereas the smart meter can recover *Message* and may want to recover the message sent by individual devices. For that purpose, it has to collude with the sink node, which is not possible under the threat model. Hence, the suggested technique protects aggregated data privacy.

### **7. Performance Evaluation**

In this section, we analyze the communication and processing overhead of household appliances, sink nodes, and smart meters.

### *7.1. Communication Overhead*

The proposed privacy preserving scheme aggregates the data from different subsets into one and the smart meter can recover the mean and variance of the individual subset. To demonstrate the efficiency of the proposed technique, we compare it to the basic Paillier encryption [80] in which the bit length of *n*2 is 2048 and that of n is 1024. Therefore, the communication overhead from N devices to sink node is 2048 × N bits, because each home appliance encrypts both *xi* and *x*2*i* into one cipher text, as shown in Equation (9). However, because the data in the BPE is encrypted into two cipher texts, the transmission cost is doubled, and the overhead is 4096 × N bits.

In the proposed scheme, the sink node and smart meter are independent of the number of devices as the aggregation is undertaken at the sink node. In BPE, the communications cost is dependent on the number of subsets. If there are k subsets, then the communications cost is 4096 × k bits, while in the proposed scheme, all data have been aggregated in one cipher text. Therefore, the communications overhead from the sink node to the SM is only 2048 bits. Figure 2 plots a graph for the communications overhead from the home appliances to the sink node, and Figure 3 shows the communications overhead pattern of the sink node to the SM. As a result, it is obvious that the suggested method is superior to BPE in terms of communication costs. We use the Chinese remainder theorem that enables a careful parameter choice and, hence, in real-time scenarios the message size is small.

**Figure 2.** Communication overhead from home appliances to sink node.


*2.2. Regulatory Overview and In-Depth Analysis of the Main Italian Legislative Measures for School Construction*

As previously specified, in Italy, since the 1970s, there have been many legislative measures in the field of school building aimed at providing planning and/or design indications for the execution of interventions on existing schools, preferring a problem-solving-based approach, and not one of integration and compensation between multiple effects deriving from the same settlement transformation intervention. The main reference standard documents currently in force (2019) for school building are the:

the management of interventions on school building at national and regional level is carried out, is established;

These standards can be used as a reference both for existing buildings and new constructions. From each one, the technical-regulatory objectives to be pursued can be defined, and the corresponding parameters (Evaluation Criteria) can be identified with which one it is possible to express the degree of achievement of each objective.

Table 1 shows some of the main regulatory references in the field of Italian school buildings. The basic purpose is described for each of them, and the corresponding evaluation criteria is specified with which one can express the level of conformity between the actual state of the school building and the legal requirements to be complied with. With respect to the list of normative references described above, in Table 1 the normative references containing provisions for Italian school field still in force are taken into consideration.

**Table 1.** Main regulatory references for Italian school buildings and corresponding evaluation criteria.

Description of the Information Contained in the Data Sheet for School Buildings

Among the normative references illustrated above, the Law n◦ 23/1996 (Norms on school building) underlines the necessity to carry out a planning phase of interventions aimed at the conservation of the

existing building considering the survey about the characteristics of the school structures in use and of the urban context. Through this law, the *Anagrafe Regionale dell'Edilizia Scolastica* (ARES) is established with the aim of systematizing the information system on the regional school property assets.

Each region and autonomous province independently manages access to ARES, and the provinces and municipalities are responsible for compiling, updating and implementing the data collection forms for each individual school building. This was done through direct surveys and inspections conducted at the school of interest.

The module for collecting data on school buildings contained the elements needed to acquire information ascertained through the completion of two questionnaires:

(i) *Questionario Edificio* (QE), aimed at collecting elements to evaluate quantitatively and qualitatively the school in use;

(ii) *Questionario istituzione scolastica*, aimed at collecting information on individual school units, i.e., whether or not there are several school units in the same building, what type and how they are organised.

The data that can be deduced from (i) on the school building are summarized in the Appendix of *Questionario Edificio* in ARES by Law 23/1996. For a more detailed specification, please refer to the Instruction Manual of QE in ARES for the compilation of the school building stock survey sheets.

In order to have a clear view of the current state of the suitability of the school building for its functions when new work is to be carried out on the building, the relevant data which are particularly interesting include the most recent work on the structure (Point 9.0—Subsequent transformations). For each intervention, it is important to establish the intervention class and the execution year. The classes of intervention referred to in the survey, aimed at determining the state of maintenance, and consequently the level of adequacy of the building to the functions for which the building is intended, both from the technological point of view and in terms of amount of space, are:

The state of maintenance regarding the building works and systems (Point 13.0—State of conservation) is evaluated qualitatively by attributing a score according to the following classification: 6 = does not require any intervention; 5 = requires partial maintenance; 4 = requires complete maintenance; 3 = requires replacement or partial refurbishment; 2 = requires replacement or complete refurbishment; 1 = requires ex-novo installation; X = system is not necessary;

The functional and dimensional characteristics of the rooms include the location, functional destination, size of the rooms on each floor of the building, shape (when the floor plan of a room is clearly different from a rectangular square, and it is difficult to carry out the educational activities), natural and artificial lighting, hygienic conditions (dependent not on poor cleaning but on a physical deficiency of the inner shell), and the acoustic conditions of each space.

Depending on the observations and specifications previously made, it is clear that in order to plan interventions on existing school buildings, it is necessary to comply with a number of regulatory requirements, and to obtain information on the building to be recovered. In general, the information concerning the school structure can be of various types. Each one regards a specific school aspect, both technical and organizational-functional, especially regarding the way in which the available space is used.

On the basis of the information system proposed by Law 23/1996, it can be observed, however, that the phase of reconnaissance of the state concerning the school building does not include the acquisition of data in the urban context.

With reference to the Constitutional Court (judgments 62/2013, 284/2016 and, lastly, 71/2018), within the discipline on school buildings "more subjects intersect, such as "territorial governance", "energy" and "civil protection". Thus, the proposed evaluation methodology, as specified in the following paragraph, aims to jointly consider regulatory aspects, school building features, and characteristics of the territory.

Some of the information illustrated in Table 2 is identified from the data list contained in *Questionario Edificio*. For each one, an alternative reference document to the QE is indicated, from which the data relating to the information to be quantified can be extrapolated. This is with a view to creating an information system useful for the planning and execution of projects that respect the characteristics of the building to be recovered, and that can also satisfy the citizens' needs in view of the distinctive characteristics of the local economy expressed in terms of existing services and those potentially settled in the area. In order to establish what further activities can be offered to the community by taking advantage of the internal spaces and connected externally to the school, it is necessary to identify the users' class to which the type of service provided is allocated (as specified in the following Section 3.1.2).

### **3. Evaluation Framework**

The proposed evaluation method seeks to verify the technical, regulatory and management conditions of existing school buildings (especially with regard to the use of available space for educational activities) in order to plan transformation and/or conservation activities, not only observing the regulatory requirements in force and the physical characteristics of the building, but also the needs of the community.

The method consists of two steps: (a) knowledge phase, (b) evaluation phase. In the first phase, the reference legislation is analyzed, and technical-management information on the school is collected. In the second phase, the regulatory requirements are verified, and the main methods of intervention are defined. This is with regard to both the physical and management system of the building and the services that the school can offer to the community.

Figure 1 shows the outline of the proposed evaluation methodology and highlights the mutual relations between each of its sub-phases. The diagram in Figure 1 aims to graphically illustrate how the physical and functional characteristics of the school buildings are checked and evaluated in accordance with the reference regulations, as well as to program interventions compatible with existing community's needs. Compatibility and consistency assessments provide useful information on the types of interventions to be implemented for the renovation/recovery of the building.

Each step (knowledge phase and evaluation phase) is made of specific sub-phases analysed in the following.


**4. Findings**

#### *4.1. Quantitative Findings*

Table 2 depicts the number of publications per journal. We selected journals with two or more articles published on non-financial information regulations because the majority (35 journals) published one study from the selected years. *Sustainability* was the journal with the highest number of publications, with seven papers. It published articles on the issue of non-financial information regulations four years in a row, from 2017 to 2020. Journals published on disclosure laws on non-financial information do not necessarily have to be pure legal journals. The journals include *Accounting in Europe*, *Business Strategy and the Environment*, *Journal of Business Ethics*, *Meditari Accountancy Research* and *Sustainability*. The results implied that the requirement of disclosures on non-financial information is also a concern and of interest to other disciplines beyond the legal fraternity. This trend indicates the importance of disclosing non-financial information and how to manage it to ensure its effectiveness. The findings also highlighted the disciplines interrelating with disclosure laws on non-financial information, including accounting, business management, ethics, and the economy.

**Table 2.** The total number of publications per journal from 2014 to early 2021.

Source: WoS and Scopus database.

Figure 2 shows that the number of papers discussing non-financial information regulations generally increased from year to year. However, only one article was found to be published in 2021, as this review paper covers publications until March 2021. Nevertheless, seeing the trend in 2020, the authors anticipated that more publications would follow on this issue as the year progressed, considering its importance.

**Figure 2.** Number of documents by year.

Findings on the publications according to country can assist other researchers in identifying the researched areas and determining future research areas [43]. Based on the findings in Figure 3, researchers focused on the issue of non-financial information regulations in China, European Union countries in groups, Italy, the United Kingdom, and the United States. The publication containing studies on the European Union countries in groups continued every year from 2014 to 2020. However, in 2017, no study related to this matter was conducted in European Union countries. Nevertheless, studies were conducted in European Union member-states individually; Italy, Poland, and Romania, in 2017. The study on European Union countries in 2014 unsurprisingly discussed the integration of sustainability information in corporate reporting. It analyzed how legislation could support integrated reporting [44], considering that the introduction of < IR > Framework was in December 2013. The study conducted in China in 2014 investigated the relationship between regulatory pressures exerted by the government on the breadth of corporate social-responsibility reports [45].

**Figure 3.** Documents by country.

There were two studies conducted in the Netherlands in the year 2014. The studies deliberated on using extraterritoriality and soft law in admonishing companies that disregarded human rights in their operation [46], and examined the effect of the regulatory environment on the application of integrated reporting [44]. The study conducted in Africa in 2014 intended to make companies accountable for false and misleading statements voluntarily published by them [47]. Interestingly, there has been no study on non-financial information regulations in Africa for four years consecutively, after 2014. The most probable reason for this absence lies in the fact that in that most studies conducted in Africa focused on examining the acceptance of integrated reporting. Examples include works that examined first, the development of integrated reporting by companies [48,49]; second, the challenges in preparing integrated reports [50]; third, the potential of integrated reporting in developing value-creation [51]; and fourth, the application of <IR> principles to sustainability disclosures [52], rather than scrutinizing the impacts of regulations on non-financial disclosures. A study conducted in the United Kingdom in the year 2014 focused on the area of corporate accountability regarding human-rights violations [46].

In 2015, the number of publications relating to non-financial information regulation increased to 11. Nevertheless, they then went down to only five publications in 2016, where the studies were conducted in Canada, China, European Union countries in groups, India, and the United States. However, the number of publications relating to this matter surged to 24 in 2019. There were four studies in the United Kingdom and three in both Italy and the United States. Therefore, inspecting the issues explored and discussed by these

studies is interesting. A study in the United Kingdom focused on corporate accountability for labour standards, precisely the Modern Slavery Act 2015 [53]. The subsequent study examined the effects of mandatory disclosure under the Companies Act 2006, Regulations 2013, on disclosed sustainability-information [54]. Another study examined how the United Kingdom adopted the 2014/95/EU Directive into their legal system [55], while Ho and Park [29] investigated the approaches by the United Kingdom regarding ESG matters.

Studies in the United States investigated the disclosure approaches to ESG matters [29], corporate accountability of companies that failed to disclose critical climate-change information [56] and the potential of disclosure mechanisms in ensuring the preservation of human rights by companies [57]. As for the studies in Italy, the researchers explored the level of the legal-system adaption to the requirements of the 2014/95/EU Directive [55], the quality of non-financial information after the implementation of the Legislative Decree 254/2016 [58] and the effects of the introduction of mandatory non-financial-reporting obligations by the 2014/95/EU Directive [27]. It is noteworthy that the studies conducted in the United Kingdom and the United States overlapped in one area relating to corporate accountability, implying concern by the scholars.

Another significant finding presented in Figure 3 is a lack of studies in countries aside from the United Kingdom, the United States, European countries, and China. Only three studies relating to disclosure laws were conducted in Southeast Asian countries. Mohan [59] supports the National Action Plans by Southeast Asian countries generally to take cognizance of the international standards for uniformity among these countries, especially on corporate accountability. A study in Singapore examined the effects of mandatory sustainability-reporting-requirements on the stock market [60]. A study in Indonesia investigated the relationship between government regulations and the extent of sustainability-information disclosure [61]. There has been no study in Brunei, Myanmar, Cambodia, Timor-Leste, Laos, Malaysia, the Philippines, Thailand, or Vietnam. Therefore, researching these countries regarding disclosure laws on non-financial information would be of significance and interest.

#### *4.2. Qualitative Findings*

The second part of the findings focuses on the qualitative discoveries derived from the ATLAS.ti version 9 software. Eight themes or patterns have been identified in trying to answer the research question. The themes identified are fiduciary duties of directors, corporate accountability, disclosure approaches, stakeholder engagement, the effectiveness of regulatory interventions, the impact of regulations, Directive 2014/95/EU, and the role of different actors, as presented in Figure 4.

Despite being subject to the same category, there are different groupings for the themes of Directive 2014/95/EU, the impacts of regulations, and regulatory interventions' effectiveness. The basis for such differentiation lies in the fact that Directive 2014/95/EU attracts a specific discussion. On the other hand, the other two themes touched on other disclosure laws, such as the United Kingdom Modern Slavery Act 2015 and the Bribery Act 2010. The articles grouped under the impact of regulations discussed the effect of regulations on different constructs, such as quality of reporting, the quantity of reporting, and corporate behaviour. The articles grouped under the theme of the effectiveness of regulatory interventions discussed the extent or effectiveness of a specific disclosure law in administering non-financial disclosures in corporate reporting. Several articles were grouped under more than one theme. The qualitative findings then proceeded as follows, for the analysis of each theme.

**Figure 4.** The patterns and trends in the literature on non-financial information regulations. Source: Authors' interpretation from the documents synthesized from WoS and Sco-pus database.


**4.3 Validation of opposing jet configuration**

Numerical results are compared with the experiment by NASA Langley in 1987 [24]. The radius of the cylinder is 25.4 mm and the thickness of the steel is 12.7 mm [35] (**Figure 9**). The free flow parameters are as follows (**Table 5**).

To investigate the effect of thermal protection systems we carried out a comparison between the opposing jet model and spike model at 26 km under the following free stream conditions illustrated in **Table 6**.

### **Figure 9.**

*Flow field wave structure of opposing jet; (a) our numerical photograph result; (b) experimental Hayashi photograph [24].*

#### **Table 5.**

*Opposing jet free stream boundary condition [35].*


*Article* **A Miniaturized Device Coupled with Digital Image Correlation for Mechanical Testing**

**Daniel J. Cruz 1 , Jose Xavier 2,3,\* , Rui L. Amaral 1 and Abel D. Santos 1,4**

**Abstract:** Miniaturized mechanical testing based on small sample testing technology is a powerful technique to characterize the mechanical properties of different materials, and it is being used in different application fields. However, the small size of the specimens poses several challenges because the results are highly sensitive to measurement accuracy and the corresponding mechanical properties can change substantially due to the so-called specimen size effect. In this work, a novel testing device based on miniaturized specimens is presented. The equipment is designed to test materials in tensile and compressive loadings, but it is also capable of performing reverse-loading tests. Buckling of the specimen is an inherent phenomenon in compression loadings, especially for thin materials. Therefore, specimen geometry is properly studied and optimized to mitigate this effect. To evaluate the deformation of the specimen, the digital image correlation (DIC) technique is used to capture the full-field strain in the central gauge section of the sample. A sensitivity analysis of the DIC setting parameters was performed for this application. To evaluate the performance of the developed system, experimental results of monotonic tests and tests with reverse loadings (tension-compression) are presented, considering two high-strength steels (DP500 and DP780).

**Keywords:** miniaturized specimen; digital image correlation (DIC); Small Specimen Test Techniques (SSTT); Bauschinger effect; sheet metal characterization

#### **1. Introduction**

The mechanical behavior of engineering materials is typically determined by experimental tests using the so-called "macro specimens", with dimensions on the scale of millimeters. However, in some applications, the material to be characterized is not available in the form of large blocks [1], which has recently motivated an increasing interest in evaluating the mechanical properties of materials with "miniaturized" specimens, also known by other names such as "small-size", "sub-size", or "mini" specimens [2,3]. Miniaturized specimens were first developed in the 1970s to characterize materials for in-service nuclear reactor pressure vessels [4]; this is an application in which surveillance capsules are used to periodically monitor the degradation of vessel materials. These capsules typically have limited storage space, and for this reason, a limited specimen size is required [5]. As a result, Small Specimen Testing Technology (SSTT) has been developed and applied in different application fields such as residual life assessment of in-service components [6], local measurement of properties in the welding industry [7,8], the additive manufacturing process [9,10], and the evaluation of mechanical properties of miniature devices such as microelectromechanical systems (MEMS) [11,12].

Several non-standard test techniques are currently available to evaluate the mechanical behavior of materials using miniaturized specimens [3]. One of the most used SSTT is

**Citation:** Cruz, D.J.; Xavier, J.; Amaral, R.L.; Santos, A.D. A Miniaturized Device Coupled with Digital Image Correlation for Mechanical Testing. *Micromachines* **2022**, *13*, 2027. https://doi.org/ 10.3390/mi13112027

Academic Editors: Luigi Sirleto and Giancarlo C. Righini

Received: 4 November 2022 Accepted: 17 November 2022 Published: 19 November 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

the miniaturized tensile test (MTT) [13], which uses specimens with dog-bone geometry with nominal dimensions ranging from several millimeters to hundreds of microns. No standards have yet been established for MTT, which has resulted in different specimen designs. The design of miniaturized specimens currently uses one of two approaches: (1) simply scaling down from standard specimens, or (2) scaling down from standard specimens, with a subsequent adjustment of some dimensions. A comprehensive review of this topic can be found in the recent work of [2].

Miniaturized tensile testing encompasses several challenges related to specimen manufacturing, handling, and procedural issues. In addition, the miniaturization may cause the so-called "specimen size effect", which leads to different behaviors at the miniaturized and macro scales. Specimen size effects can be categorized into two different levels. At the first level, the test results obtained from miniaturized specimens reflect the bulk mechanical properties described by continuum mechanics. However, for some materials, the results may deviate from standard specimens, and additional data processing techniques, e.g., the inverse Finite Element Method (FEM), are required [14–16]. On the other hand, at the second level, the specimen can no longer reflects the bulk properties as the thickness or diameter is further decreased. Typically, the critical thickness value is defined by the ratio *t*/*d* between thickness, *t*, and grain size, *d*, of the material. Generally, miniaturized specimens with a *t*/*d* ratio between 6 and 10 can achieve standard bulk behavior.

Sheet metal forming is an important manufacturing process in the production of metal components, applied in a wide variety of industrial areas such as automotive, naval, electrical, aeronautics, and aerospace [17]. Sheet metal forming includes simple processes such as bending, stretch forming and spinning, as well as more complex ones such as roll forming and deep drawing. In these processes, material behavior assessment is crucial for the correct calibration of material constitutive models. Therefore, a wide variety of tests are established in the literature to obtain reliable experimental data used in the identification of constitutive parameters.

The mechanical characterization of metallic sheets has standardized test procedures, by using the tensile test (uniaxial loading), according to ISO 6892-1 [18] and the hydraulic bulge test (biaxial loading), according to ISO 16808 [19]. On the other hand, no standardized test setup has yet been established to characterize the compressive loading of sheet metal materials [20]. Moreover, an uniaxial compressive stress state may cause buckling in flat sheet metal materials for small compressive strains. Therefore, the test method becomes invalid due to the change on the stress state. Over the past few years, different approaches have been developed to minimize buckling during compression tests. These approaches can be categorized into three groups: (i) the pack method, which involves packing a group of samples together in the thickness direction; (ii) a single sample with a lateral support system; and (iii) miniaturized test samples with a low gauge length/thickness ratio.

The pack method, developed by Aitchison et al. [21] and later improved by Jackman [22], involves joining together a group of samples along the thickness direction to form a unique specimen. With this strategy, the length/thickness ratio is reduced, and therefore the buckling instability under compression is minimized. A variation in the pack method was used by Yoshida et al. [23], in which five specimens were glued together and then covered with clamping plates. Because the pack method requires a large number of specimens and complex preparation, some authors prefer to use a single sheet with a lateral support system. Boger et al. [24] sandwiched a single flat specimen between two solid support plates using a constant restraining side force. To prevent buckling in the thickness direction (T-buckling), width direction (W-buckling), and unsupported gap (L-buckling), an oversized dog-bone specimen was designed. A device with two sets of fork-shaped supports was proposed in the same year by Sekine and Kuwabara [25] to eliminate the uncovered area in tensile-compressive relative movements. Because the male and female dies slide past one another as the sample is compressed, this design made it possible to support the entire length of the specimen.

The gauge-length/thickness ratio of the specimen has the biggest impact on the buckling phenomenon in compression. In order to effectively study the compression behavior of different sheet metal materials, test methods using miniaturized specimens have recently been suggested. Tritschler et al. [26], tested a titanium alloy with a thickness of 1.2 mm in compression using specimens with dimensions of 2.5 mm in length by 2 mm in width, achieving about 2% compressive strains. Hußnätter [27] used a specimen geometry of 2 mm length by 2 mm width to characterize the compressive behavior of magnesium alloy AZ31, achieving about 5% compressive strains before the specimen buckles.

Cruz et al. [28] developed an experimental equipment to perform uniaxial mechanical tests on miniaturized metallic specimens. The current work is a step forward in these developments, exploring the applicability of this equipment to describe the mechanical behavior of sheet metal materials, not only for monotonic uniaxial tests, i.e., tension and compression, but also for reverse loading tests, i.e., tension followed by compression or vice versa. Moreover, the digital image correlation technique was coupled with the miniaturized test device in order to quantify the quality of the test method regarding the specimen deformation at the scale of observation, in a balance between spatial resolution and accuracy [29].

#### **2. Materials and Methods**

#### *2.1. Materials*

In the present work, advanced high-strength steels (AHSS) produced by SSAB company were used, namely two grades of dual-phase (DP) steel sheets (DP500 and DP780) with an initial thickness of 0.8 mm. The dual-phase (DP) steel is an advanced high-strength steel (AHSS), with a good balance of strength and ductility due to its martensite islands in a ferrite matrix, as represented in Figure 1. Table 1 lists the percentage of alloying present in each material chemical composition.

#### *2.2. Design of an Uniaxial Test Equipment for Miniaturized Specimens*

The experimental test procedure for miniaturized specimens followed a similar procedure to the one typically carried out for conventional standard specimens, which is based on moving grips and results of force as a function of elongation/strain data. However, a very accurate procedure is required because results can be significantly affected by the testing protocol, including sample preparation, test setup, experimental operation, measurements, and data processing [2]. Several authors have emphasized different aspects

to be taken into account when testing miniaturized specimens [3]: (a) the importance of the specimen-fabrication technique, which can have consequences for the measured properties, particularly on the yield strength; (b) good control and accuracy of the specimen dimensions, which directly relate to the reliability and reproducibility of the test data; (c) high-quality measuring instruments; and (d) the importance of a precise and well-aligned test frame, in order to minimize the extraneous bending strains.

Based on these considerations, a new type of equipment was developed, specifically designed to test miniaturized specimens on a reduced scale, and we refer to it as a Miniaturized Specimen Tester Device (MSTD) [28]. As represented in Figure 2, MSTD consists of three groups of components with their well-defined function. Firstly, there is the "motor group", which includes a ball screw, driven by a stepper motor, which guarantees continuous and precise movement of the specimen, making it possible to guarantee a minimum displacement of 2 µm for each step of the motor. Secondly, there is the "frame group", which ensures the layout and guidance of the various components, consisting of four linear bearings that slide on two ground shafts. Finally, there is the instrumentation group, consisting of a 5 kN "S"-type load cell and two extension measurement components: a linear ruler for measuring and controlling the positioning of the moving table and a digital image correlation (DIC) system. The MSTD's layout allows all components to be on a common plane, which is advantageous for a precise alignment of the specimen during the test. Consequently, the appearance of undesirable bending moments is minimized, especially in compression, which is also essential to prevent the buckling of the specimen. The equipment, designed to test miniaturized specimens, is capable of withstanding a maximum axial force of 2.5 kN.

The preparation of the specimen is crucial in MTT because miniaturized specimens are more susceptible to large-dimensional deviations. Accordingly, in this work the specimen profile was obtained by wire electrical discharge machining (EDM) using a 0.25 mm wire diameter. This process is widely used for miniaturized specimen manufacturing, resulting in good dimensional accuracy and surface finishing [2].

**Figure 2.** Miniaturized Specimen Tester Device (MSTD): (**a**) setup, and; (**b**) geometry of the miniaturized specimen.

#### *2.3. Testing Conditions*

In order to evaluate the monotonic behavior of DP500 and DP780 steels, tests were performed with tensile and compressive loading using miniaturized samples and the geometry shown in Figure 2b. Additionally, the same materials were tested on macrotensile test pieces, in accordance with ISO 6892-1 [18]. These macro standard tensile tests were performed on a 300 kN Instron 5900 R testing machine. The experimental conditions for both geometries are shown in Table 2. The specimens were tested at a constant test speed, resulting in an initial strain rate of approximately 10−3 s −1 , until the material failed. To ensure the repeatability of the results, two samples were systematically tested for each geometry.

**Table 2.** Experimental conditions for macro and miniaturized samples.

#### *2.4. Digital Image Correlation*

#### 2.4.1. Choice of the Optical Technique

The DIC technique was coupled with the mechanical tests to quantify the deformation of the material across the central region of interest for both macro and miniaturized specimens, with a suitable spatial resolution and accuracy at the scale of observation. Traditionally, in carrying out experimental tests (e.g., macro tensile test [18]), the strain components are measured by means of a strain gauge attached to the sample. However, for miniaturized specimens, these instrumentation devices may not be suitable due to the reduced sample size. Nevertheless, a non-contact optical technique can be conveniently used to access full-field strains over the target region. For that purpose, several white-light and interferometric optical methods have been proposed in the literature [30]. Among them, DIC has been increasingly used in recent years in the field of experimental solid mechanics [31,32]. This technique was introduced in the 1980s [33,34] and allows the measurement of full-field kinematic quantities by tracking unique features on images taken at different deformation stages, typically by means of a subset-based approach [35–38]. The DIC technique uses a camera-lens optical setup, and it allows variable sensitivity and resolution in strain measurements. Since the DIC method uses digital images, several high spatial-resolution digital devices can be used for image acquisition. For example, DIC can be coupled with Optical Microscopy (OM) [39,40], Scanning Electron Microscopy (SEM) [41,42], and Scanning Tunneling Microscope (STM) [43] to realize microscale-tonanoscale deformation measurements. In the 2D DIC set-up, the planar surface of a specimen is imaged stringently parallel to the camera sensor and the surface must have a random pattern of a contrasted speckle pattern. The creations of a quality pattern for DIC at reduced scales can be challenging [44], and the accuracy of the measurements depends significantly on the spatial resolution of the imaging system [35].

#### 2.4.2. Optical System and Speckle Pattern

A 2D DIC set-up using a 5 MP CMOS digital camera was selected in this work. The mounted lens was selected with regard to the specimen size. The details of the camera-lens system used for the miniaturized tests are summarized in Table 3. In this case, conversely, to a standard 50 mm focal lens (Fujinon HF50HB-1B 50 mm lens, f/2.3 C–Mount) used in macroscopic tensile tests, a telecentric lens was selected. It had the advantage of keeping the image magnification unchanged during the mechanical test, avoiding artifacts due to, for instance, Poisson effects or slightly out-of-plane movements during the test. In this case, the undeformed image was focused by fixing the working distance at 63.3 mm with a magnification factor of 1×. Images were recorded during the tests at an acquisition frequency of 5 Hz. The central part of the specimens was initially painted by spraying a thin layer of white paint in order to make the diffuse reflectivity of the surface uniform. The final speckle pattern was then created by spreading dark paint using an airbrush (Harder & Steenbeck—Evolution Silverline Solo) with a nozzle set of 0.2 mm. Considering the region of interest and camera-lens system, an average speckle size of about 6 pixels (21 µm) was obtained in the speckle quality evaluation. According to international guidelines [45], this means that a subset size larger than 21 pixels (0.075 mm) or 31 pixels (0.111 mm) must be used in the DIC measurements for accuracy.

**Table 3.** Camera-lens optical system for the miniaturized tests.

#### 2.4.3. Selecting DIC Setting Parameters

In performing DIC measurements, a convergence analysis must be carried out for the inherent parameters in the numerical imaging method [29,46]. The set of DIC parameters plays an important role in the spatial resolution and accuracy associated with both displacement and strain measurements [47]. A parametric analysis was therefore carried out in this study, using the Performance Analysis Tool of the MatchID subset-based 2D DIC software [48]. The parametric tool allows one to carry out several DIC analyses over the same set of images by considering several combinations of parameters, as reported in Table 4. Each relevant parameter was considered over a suitable design space by defining the minimum and maximum values and a step of an integer increment (i.e., a subset size with a minimum of 21 × 21 pixels2 , a maximum of 91 × 91 pixels2 , and increment of 10 pixels, Table 4). Results in terms of signal evaluation at the center of the specimen can then be discussed in light of the preselected DIC setting parameters. The linear strain component in the *x* direction (*εxx*) was used as a reference, since uniaxial tests were carried out along this axis. Each set of DIC parameters can be associated with a corresponding spatial resolution. To quantify this metric, the following Virtual Strain Gauge (VSG) measure was used [45,48]:

$$\text{VSG} = [(\text{SW} - 1) \times \text{ST}] + \text{SS} \quad \text{[pixels]} \tag{1}$$

where SW stands for a strain window, i.e., the number of displacement data points used to fit a local polynomial in the neighbor of the subset centroid under evaluation; ST is the subset step; and SS is the subset size. The VSG can be expressed in physical units of millimeters if the image conversion factor of the optical system is taken into account (Table 3).

In order to evaluate the *εxx* strain signal with regard to the DIC settings as a function of the VSG, three points of reference at the centre of the specimens were systematically analyzed for both tensile and compression tests, for the DP500 steel, as schematically shown in Figure 3.

Figure 4 shows the evolution of *εxx* as function of points A, B, and C (Figure 3) for both tensile and compression tests at two different deformation stages corresponding to low and high average strains, respectively. As can be observed, the DIC setting parameters converge systematically to an average reconstruction of the linear strain. In a more local analysis, Figure 5 shows the results of the reconstructed *εxx* value at point B in the center of the specimen (see Figure 3) as a function of VSG, for the different DIC setting parameters for both tensile and compression tests and for two stages at low and high strain values, respectively. A good balance is to be found with regard to the spatial resolution in reconstructing the strain field. As can be seen in all scenarios, the signal is rather stable with regard to a set of parameters. This may be expected for a homogeneous material and uniaxial loading configuration.

**Table 4.** Parameters used in the Performance Analysis Tool of MatchID 2D DIC.

\* Local least-squares fitting approach for the strain evaluation.

(**a**) (**b**) **Figure 3.** Points of reference used in evaluating the reconstruction of the *εxx* strain signal with regard to the DIC settings and as a function of the VSG for both (**a**) tension and (**b**) compression tests.

Figure 6 shows a complementary analysis plotting the *εxx* value along the coordinate *x* at the center of both the tensile and the compression tests for several VSG values. As can be seen, the DIC setting corresponding to higher values of VSGshows a smooth signal evaluation, as it may be expected, in contrast with lower VSG DIC settings. In conclusion, from the analysis in this work, the DIC settings reported in Table 5 were finally selected in a suitable balance between spatial resolution and accuracy of results for the current application.

**Table 5.** DIC and strain settings selected in the full-field measurements.

**Figure 4.** Evaluation of *εxx* for points A, B, and C (Figure 3) as a function of DIC parameters expressed with regard to Virtual Strain Gauge (VSG) for both tension and compression tests, at two stages of deformation corresponding low and high strains. (**a**) Point A—Tension; (**b**) Point A—Compression; (**c**) Point B—Tension; (**d**) Point B—Compression; (**e**) Point C—Tension; (**f**) Point C—Compression.

**Figure 5.** Evaluation of *εxx* at point B (Figure 3) as a function of DIC parameters expressed with regard to Virtual Strain Gauge (VSG). (**a**) Tension-Low Strain; (**b**) Tension—High Strain; (**c**) Compression— Low Strain; (**d**) Compression—High Strain.

**Figure 6.** Evaluation of *εxx* as a function of the coordinate *x* for several Virtual Strain Gauge (VVSG) values for both tension and compression tests.

#### **3. Results and Discussion**

*3.1. Monotonic Tension and Compression Tests*

For the reconstruction of the strain signal during the mechanical tests, an average value was systematically calculated at the centre of the specimen corresponding to a region of 2 mm × 2 mm. Figure 7 shows a representative stress–strain curve for both DP500 and DP780 steel materials. The engineering stress–strain curves, shown in Figure 7a, indicate that the results obtained from both macro and miniaturized specimens are quite

comparable. In fact, for both materials, the work hardening is very similar for the mini and macro geometry, as seen in the true stress–true strain curve (*σ* − *ε*) shown in Figure 7b. However, it can be seen that after the point of ultimate tensile strength, the miniaturized geometry presents a superior total elongation when compared to the macro geometry. This difference can be justified by the local necking that occurs after reaching the point of maximum stress, which translates into a localized strain zone whose size is independent of the gauge length. Therefore, for smaller gauge lengths, the localized strain has a higher impact on the total elongation value. Additionally, the variation in necking behavior should be taken into account as well. As reported by [49], total elongation and post-necking elongation tend to increase with a decreasing ratio, *w*/*t*, between specimen width and thickness. In fact, when *w*/*t* decreases, the shear deformation may occur in both the upper (lower) and side planes, resulting in complex stress distribution, and the localized necking changes to diffuse necking, which leads to higher post-necking elongation.

**Figure 7.** Experimental tensile stress/strain behavior for DP500 and DP780 steels: (**a**) engineering stress–strain curves; (**b**) true stress–strain curves.

The hardening behavior of the two dual-phase steels under compressive loading was also evaluated using the miniaturized specimens. The true stress–true strain curve (*σ* − *ε*) for each material is shown in Figure 8a. Additionally, this figure represents, for each material, the tensile hardening curve obtained using miniaturized samples (curve defined by points). For an easy comparison with these two tensile curves, the compression results were reflected for the first quadrant and represented as a dashed line. When analyzing results, it was verified that the work hardening was similar in the tension and compression for both materials; however, in terms of compression, the maximum strain values in the compression test were lower than those reached in tension. As seen in Figure 8b, for a strain value close to 7%, the force reaches a maximum value. At this point, the test becomes non-planar and the buckling phenomenon becomes evident. This instability is translated in the decay of the stress, which means that the test is not longer valid.

#### *3.2. Reverse Loading Tests: Tension–Compression*

As represented in Figure 9, a typical hardening curve of a metallic material during reverse loading can be divided into four interconnected phenomena [50]: (a) the Bauschinger effect, which describes a reduction in the yield stress after the occurrence of plastic deformation during the initial loading; (b) the transient behavior, which is given by a smooth elasto-plastic transition characterized by a change in the work-hardening rate; (c) the work hardening stagnation, given by a transient behavior in which a stress–strain plateau is observed; (d) permanent softening, characterized by lower levels of stress after the transient period when compared to monotonic loading at a given accumulated plastic strain. These phenomena will hereafter be referred to as Bauschinger behavior. The experimental characterization of Bauschinger behavior is crucial to calibrate advanced constitutive models, e.g., the multi-surface Yoshida–Uemori (Y–U) [51] kinematic hardening model and the

Homogeneous anisotropic hardening (HAH) model proposed by [52], and consequently improves the accuracy of numerical simulation results of sheet-metal-forming processes. Different authors have shown that the Bauschinger behaviour has a strong influence on the prediction of the springback [53–55]. This phenomenon, which refers to the elastic distortion of a part when forming tools are removed, is a major issue in sheet metal forming processes and its mitigation is crucial.

**Figure 8.** True stress–strain curves from monotonic compression tests and comparison with tension tests; (**b**) buckling point at the compression test.

**Figure 9.** Schematic illustration of stress–strain response during forward and reverse deformation.

In this work, the Bauschinger behavior of the DP500 and DP780 materials was evaluated through uniaxial tests with reverse loading (tension–compression) using miniaturized specimens. The test includes a step of tensile loading, followed by a compression step. Initially, the tensile step is performed up to a defined pre-strain, *ε* ∗ *p* . In this sense, three pre-strain values were selected, 0.025, 0.050, and 0.075. After this pre-strain in tensile loading, the compressive step was carried out until the buckling phenomenon began to be evident, which is shown, as already mentioned, by a decay of the stress value in the stress–strain curve. The results for DP500 and DP780 steels are shown in Figure 10a,b, respectively. For each material, the monotonic uniaxial tensile curve is also represented for comparison purposes. It can be seen that the tension step of tension–compression test

matches the uniaxial tensile curve very well, for both materials, which confirms that the equipment provide excellent repeatability of results. The compression step curves for the three different pre-strain values show a similar changing trend across the entire strain range, indicating that the device can provide reliable results for the studied materials. To analyze the Bauschinger behavior of the dual-phase steels used, the compression-step-hardening curves were reflected to the first quadrant, as shown in the dashed line in Figure 10. The yield point for the uniaxial tensile curve and the reflected reverse loading curves, represented in round points, were taken at the 0.2% plastic strain and are summarized in Figure 11a. It can be observed that the studied materials present a strong Bauschinger effect, since the reverse yield stress is significantly reduced during reverse loading. A strong permanent softening effect is also visible, since the curves do not coincide with the monotonic tensile curve and tend to be approximately parallel.

**Figure 10.** True Stress–Strain curves from tension–compression tests for (**a**) DP500 and (**b**) DP780.

**Figure 11.** (**a**) Yield stress, *σY*, and (**b**) Bauschinger Ratio, *β*, values for different pre-strains, *ε* ∗ *p* .

To quantify the Bauschinger effect, several indicators have been proposed over the years. One of the most used, the so-called Bauschinger ratio, *β*, is based on stress values and is expressed mathematically by Equation (2):

$$\beta = \frac{\sigma\_{\text{max}}^T + \sigma\_{Yield}^C}{\sigma\_{\text{max}}^T} \tag{2}$$

where *σ T max* represents the maximum stress in the forward load (tensile step) and *σ C Yield* represents the yield stress in the backwards or reverse load (compression step)-point B and C, respectively, in Figure 9. The Bauschinger ratio quantifies the magnitude of the Bauschinger effect. The higher the parameter, the greater the Bauschinger effect of the material. In the case of isotropic hardening, the Bauschinger ratio is equal to zero. The Bauschinger ratio for the materials under this study is represented in Figure 11b, which is dependent on the applied pre-strain. From these results, it is possible to observe that both materials have a strong Bauschinger effect, since whenever the loading direction is reversed, the yield stress is reduced by about 70–80% compared to the value of the stress in the opposite direction.

#### **4. Conclusions**

This work presents the development of equipment dedicated to performing uniaxial, tensile–compression tests on miniaturized metallic specimens. The experimental equipment, called MSTD (Miniaturized Specimen Tester Device), allows the testing of specimens on a miniaturized scale, which opens the field of applications, one of them being compressive loading tests, having the advantage of overcoming buckling effects. In addition to uniaxial tensile, compression, and cyclic tests, it is also possible and expected that equipment will adapt to other fundamental tests for metallic sheets, such as shear and ductile damage tests.

The digital image correlation (DIC) technique was used to capture the full-field strain in the gauge section of the sample. A sensitivity analysis of the DIC setting parameters permitted to quantify the influence on results for subset size, subset step, strain window, shape function, and polynomial order. These three first parameters are included in the Virtual Strain Gage variable (VSG) and this variable was tested for different scenarios in tensile and compressive loading, giving a rather stable signal and a fast convergence of results. A good balance was found using a VSG of 151 px, which corresponds to 0.53 mm, thus giving the best ratio for efficiency and accuracy of results.

In order to validate the developed solution, several experimental tests were carried out using two dual-phase steels: DP500 and DP780. The tests performed on the developed equipment allowed the comparison of results for different types of loadings. Initially, monotonic tensile tests were carried out, and the results obtained with the miniaturized specimens were shown to be quite similar to those obtained by using the standardized geometry of tensile tests at the macro scale. The results of the compression tests showed that it is possible, with the developed equipment, to characterize materials up to 7.5% strain in compressive loading, without buckling, thus creating the possibility of sheet metal characterization for stress differential and Bauschinger behavior. Additionally, it was seen that the work-hardening of tested dual-phase steels (DP500, DP780) in compressive loading is very similar to their work hardening in tensile loading. Finally, uniaxial tests with reverse loadings showed that both DP500 and DP780 present a strong Bauschinger effect, since the yield stress undergoes a considerable decrease with reverse loading.

**Author Contributions:** Conceptualization, D.J.C., R.L.A. and A.D.S.; methodology, D.J.C., J.X., R.L.A. and A.D.S.; software, D.J.C., J.X., R.L.A. and A.D.S.; validation, D.J.C., J.X., R.L.A. and A.D.S.; formal analysis, D.J.C., J.X., R.L.A. and A.D.S.; investigation, D.J.C., J.X., R.L.A. and A.D.S.; resources, R.L.A. and A.D.S.; data curation, D.J.C., J.X., R.L.A. and A.D.S.; writing—original draft preparation, D.J.C., J.X., R.L.A. and A.D.S.; writing—review and editing, D.J.C., J.X., R.L.A. and A.D.S.; visualization, D.J.C., J.X., R.L.A. and A.D.S.; supervision, R.L.A. and A.D.S.; project administration, R.L.A. and A.D.S.; funding acquisition, R.L.A. and A.D.S. All authors have read and agreed to the published version of the manuscript.

**Funding:** Authors gratefully acknowledge the funding of Project POCI-01-0247-FEDER-069933- AMOB CUT3D, POCI-01-0247-FEDER-068492-AARM 4.0, POCI-01-0247-FEDER-047038-EcoForG4.0, cofinanced by Programa Operacional Competitividade e Internacionalização (COMPETE 2020), through Fundo Europeu de Desenvolvimento Regional (FEDER). The first author is also grateful to the FCT (Portuguese Foundation for Science and Technology) for the Doctoral grant 2022.12657.BD under the program POCH, cofinanced by the European Social Fund (FSE) and Portuguese National Funds from MCTES. The second author is also grateful to FCT by the support under the project UIDB/00667/2020 (UNIDEMI).

**Data Availability Statement:** Not applicable.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **Abbreviations**

The following abbreviations are used in this manuscript:

#### **References**


*2.3. EOMs: The Discrete Index-3 Form*

To transform the second order differential equation into first-order differential equations, it is common to introduce the velocity variable *q*˙ = *v*, allowing to write Equation (4) in a first-order form as:

$$\begin{cases} \mathcal{g}(\dot{\boldsymbol{\upsilon}}, \boldsymbol{\upsilon}, \dot{\boldsymbol{\eta}}, q, \lambda, \boldsymbol{\mu}) = \begin{cases} \mathcal{g}\_1 = \boldsymbol{\upsilon} - \dot{\boldsymbol{\eta}} = \boldsymbol{0}\_{\boldsymbol{\upsilon}} \\ \mathcal{g}\_2 = \mathcal{M}(q)\dot{\boldsymbol{\upsilon}} + f\_{\text{nl}}(\boldsymbol{\upsilon}, q, \boldsymbol{\mu}) + \boldsymbol{J}^T \lambda = \mathbf{0}\_{\boldsymbol{\eta}} \\ \mathcal{g}\_3 = \boldsymbol{\phi}(q) = \mathbf{0}\_{\boldsymbol{\lambda}} \end{cases} \end{cases} \tag{8}$$

These equations represent a system of constrained-DAEs of index 3 [21]. Numerical differentiation is generally employed [22] to convert Equation (8) into a discrete form.

In this work, a first order, constant time-step Backward Differentiation Formula (BDF-1), also known as Backward Euler, is employed. This choice does not imply that other differentiation schemes cannot be applied to discretize the EOMs in time. However, the time-discretization must be consistent with the defined estimation sampling and particular attention should be paid to choosing the differentiation schemes (e.g., forward or backward) because it influences the achievement of the discrete-time E-ODE form of the EOMs, as will be discussed in Section 3.

The single step method BDF-1 can be written as

$$\begin{cases} \dot{v}\_{k+1} = \frac{1}{\hbar}(v\_{k+1} - v\_k) \\ v\_{k+1} = \frac{1}{\hbar}(q\_{k+1} - q\_k) \end{cases} \tag{9}$$

where *h* represents the constant time step size and *k* ∈ Z the *k* − *th* time instance. By substituting Equation (9) into Equation (8) the discrete-time EOMs *gd* are obtained:

$$\begin{cases} \mathcal{g}\_{d1} = v\_{k+1} - \frac{1}{h}(q\_{k+1} - q\_k) = 0\_v\\ \mathcal{g}\_{d2} = M(q\_{k+1})\frac{v\_{k+1} - v\_k}{h} + f\_{nl}(q\_{k+1}, q\_k, u\_{k+1}) + f^T(q\_{k+1})\lambda\_{k+1} = 0\_\emptyset \quad . \end{cases} \tag{10}$$

This can be summarized as

$$g\_d(v\_{k+1'}, q\_{k+1'}, v\_{k'}, q\_{k'}, \lambda\_{k+1'}, \mu\_{k+1}) = 0. \tag{11}$$

Assuming the generalized coordinates *qk* and velocities *vk* at the time instance *k* to be known and given the input *uk*+1 at the time instance *k* + 1, Equation (10) is solved for *qk*+1 and *λk*+1 by substituting *gd*1 in *gd*2, and applying a Newton-Raphson-based iterative algorithm with iteration gradient *JNR*:

$$J\_{NR} = \begin{bmatrix} \mathbb{K}\_t + \beta \mathbb{C}\_t + \gamma \mathcal{M}\_{\mathbb{H}} & J^T \\ J & \mathbf{0}\_{\lambda, \lambda} \end{bmatrix}. \tag{12}$$

Here, *Kt*, *Ct* and *Mt* are the tangent stiffness, damping and mass matrices obtained from the partial derivatives of the continuous *g*2 equations in Equation (8) evaluated at time step *k* + 1:

$$K\_{\mathbf{f}} = \left. \frac{\partial \mathcal{G}\_2}{\partial q} \right|\_{k+1}; \qquad \qquad \mathsf{C}\_{\mathbf{f}} = \left. \frac{\partial \mathcal{G}\_2}{\partial v} \right|\_{k+1}; \qquad \qquad \qquad M\_{\mathbf{f}} = \left. \frac{\partial \mathcal{G}\_2}{\partial v} \right|\_{k+1}, \tag{13}$$

and *β* and *γ* are matrix coefficients function of the defined integration rule, which are given for the BDF-1 scheme by:

$$\beta = \frac{\partial v\_{k+1}}{\partial q\_{k+1}} = \frac{\partial \psi\_{k+1}}{\partial v\_{k+1}} = \frac{1}{h} I\_q; \qquad \qquad \gamma = \frac{\partial \psi\_{k+1}}{\partial q\_{k+1}} = \frac{1}{h^2} I\_q; \tag{14}$$

$$\frac{\partial \upsilon\_{k+1}}{\partial q\_k} = \frac{\partial \upsilon\_{k+1}}{\partial v\_k} = -\beta; \qquad \qquad \qquad \frac{\partial \upsilon\_{k+1}}{\partial q\_k} = -\gamma. \tag{15}$$

This time integration scheme will be exploited in the following section to set up a suitable solver structure for use in an ODE-base estimation framework.

#### **3. An Explicit Linearized Approximation for Use of the Multibody Model in State-Estimation**

The aim of this work is to combine MB models with a Kalman filter-based estimator to concurrently estimate the system states and unknown forces of a mechanism. These presented estimators, as will be discussed in Section 4, require the linearized time-discrete model matrices. In this section, a new approach to linearize the EOMs of a MB model starting from an ID-DAE form is presented. Subsequently, the linearized system matrices are analytically assembled.

By introducing the state vector *x* ∈ R*nx* for a time step *k*

$$\mathbf{x}\_k = \begin{bmatrix} v\_k & q\_k & \lambda\_k \end{bmatrix}^T \tag{16}$$

the ID-DAE form of the EOMs in Equation (11) can be re-written as:

$$g\_d(\mathbf{x}\_{k+1}, \mathbf{x}\_k, \boldsymbol{\mu}\_{k+1}) = \mathbf{0}\_{\mathbf{x}}.\tag{17}$$

In this work we assume the function *gd* of Equation (17) to be continuously differentiable in all its variables. Therefore, an explicit discrete function *fd* locally exists by

applying the implicit function theorem. Through a first order Taylor expansion, *gd* can be approximated as

$$\mathbf{g}\_d^0 + \frac{\partial \mathbf{g}\_d}{\partial \mathbf{x}\_{k+1}} \Big|\_{0} \left( \mathbf{x}\_{k+1} - \mathbf{x}\_{k+1}^0 \right) + \frac{\partial \mathbf{g}\_d}{\partial \mathbf{x}\_k} \Big|\_{0} \left( \mathbf{x}\_k - \mathbf{x}\_k^0 \right) + \dots$$

$$\frac{\partial \mathbf{g}\_d}{\partial \mathbf{u}\_{k+1}} \Big|\_{0} \left( \mathbf{u}\_{k+1} - \mathbf{u}\_{k+1}^0 \right) + \mathcal{O}(\mathbf{x}\_{k+1}, \mathbf{x}\_k, \mathbf{u}\_{k+1}) = \mathbf{0}\_{\mathbf{x}};\tag{18}$$

where (*x*0 *k*+1,*x*0 *k* ,*u*0 *k*+1) represents the linearization set point while *g*0 *d* = *gd*(*x*0 *k*+1,*x*0 *k* ,*u*0 *k*+1) for convenience of notation. Manipulating Equation (18) and neglecting the higher order terms, it can be made explicit as

$$\mathbf{x}\_{k+1} = \mathbf{x}\_{k+1}^{0} - \left[\frac{\partial \mathbf{g}\_{d}}{\partial \mathbf{x}\_{k+1}}\Big|\_{0}\right]^{-1} \left(\mathbf{g}\_{d}^{0} + \frac{\partial \mathbf{g}\_{d}}{\partial \mathbf{x}\_{k}}\Big|\_{0} (\mathbf{x}\_{k} - \mathbf{x}\_{k}^{0}) + \frac{\partial \mathbf{g}\_{d}}{\partial \mathbf{u}\_{k+1}}\Big|\_{0} (\boldsymbol{u}\_{k+1} - \boldsymbol{u}\_{k+1}^{0})\right). \tag{19}$$

In compact form Equation (19) becomes:

$$\mathbf{x}\_{k+1} = f\_d(\mathbf{x}\_k, \boldsymbol{\mu}\_{k+1}) = \mathbf{x}\_{k+1}^0 + f\_d^0(\mathbf{x}\_k, \boldsymbol{\mu}\_{k+1}) + A\_{k+1}^0(\mathbf{x}\_k - \mathbf{x}\_k^0) + B\_{k+1}^0(\boldsymbol{\mu}\_{k+1} - \boldsymbol{\mu}\_{k+1}^0). \tag{20}$$

*A*0 *k*+1 = *∂ fd ∂x* 0 *k*+1 and *B*0 *k*+1 = *∂ fd ∂u* 0 *k*+1 are the linearized system and input matrices around the linearization set point.

Starting from Equation (19) and differentiating Equation (10), the Jacobians for the backward Euler implicit time-integrator can be computed as:

$$G\_{\mathbf{X}k+1} = \frac{\partial \mathcal{G}\_d}{\partial \mathbf{x}\_{k+1}} = \begin{bmatrix} \frac{\partial \mathcal{G}\_{d1}}{\partial v\_{k+1}} & \frac{\partial \mathcal{G}\_{d1}}{\partial q\_{k+1}} & \frac{\partial \mathcal{G}\_{d1}}{\partial \lambda\_{k+1}} \\ \frac{\partial \mathcal{G}\_{d1}}{\partial v\_{k+1}} & \frac{\partial \mathcal{G}\_{d2}}{\partial q\_{k+1}} & \frac{\partial \mathcal{G}\_{d2}}{\partial \lambda\_{k+1}} \\ \frac{\partial \mathcal{G}\_{d3}}{\partial v\_{k+1}} & \frac{\partial \mathcal{G}\_{d3}}{\partial q\_{k+1}} & \frac{\partial \mathcal{G}\_{d3}}{\partial \lambda\_{k+1}} \end{bmatrix} = \begin{bmatrix} I\_{\upsilon} & -\beta & 0\_{\upsilon,\lambda} \\ 0\_{q,\upsilon} & \gamma M\_{t} + \beta \mathcal{C}\_{t} + K\_{t} & J^{\top} \\ 0\_{\lambda,\upsilon} & J & 0\_{\lambda,\lambda} \end{bmatrix};\tag{21}$$

$$\mathbf{G}\_{\mathbf{x}\_{k}} = \frac{\partial \mathbf{g}\_{d}}{\partial \mathbf{x}\_{k}} = \begin{bmatrix} \frac{\partial \mathbf{\mathcal{X}}\_{d1}}{\partial \mathbf{v}\_{k}} & \frac{\partial \mathbf{\mathcal{X}}\_{d1}}{\partial \mathbf{q}\_{k}} & \frac{\partial \mathbf{\mathcal{X}}\_{d1}}{\partial \lambda\_{k}} \\ \frac{\partial \mathbf{\mathcal{X}}\_{d2}}{\partial \mathbf{v}\_{k}} & \frac{\partial \mathbf{\mathcal{X}}\_{d2}}{\partial \mathbf{q}\_{k}} & \frac{\partial \mathbf{\mathcal{X}}\_{d2}}{\partial \lambda\_{k}} \\ \frac{\partial \mathbf{\mathcal{X}}\_{d3}}{\partial \mathbf{v}\_{k}} & \frac{\partial \mathbf{\mathcal{X}}\_{d3}}{\partial \mathbf{q}\_{k}} & \frac{\partial \mathbf{\mathcal{X}}\_{d3}}{\partial \lambda\_{k}} \end{bmatrix} = \begin{bmatrix} 0\_{\mathbf{v},\mathbf{v}} & \beta & 0\_{\mathbf{v},\lambda} \\ -\beta \mathbf{M}\_{t} & -\gamma \mathbf{M}\_{t} - \beta \mathbf{C}\_{t} & 0\_{\mathbf{q},\lambda} \\ 0\_{\lambda,\mathbf{v}} & 0\_{\lambda,\mathbf{q}} & 0\_{\lambda,\lambda} \end{bmatrix};\tag{22}$$

$$G\_{\mathfrak{U}\_{k+1}} = \frac{\partial \mathcal{g}\_d}{\partial u\_{k+1}} = \begin{bmatrix} \frac{\partial \mathcal{g}\_{d1}}{\partial u\_{k+1}}\\ \frac{\partial \mathcal{g}\_{d2}}{\partial u\_{k+1}}\\ \frac{\partial \mathcal{g}\_{d3}}{\partial u\_{k+1}} \end{bmatrix} = \begin{bmatrix} \mathcal{O}\_{\upsilon,u} \\ \mathcal{U}\_{\mathfrak{f}} \\ \mathcal{O}\_{\lambda,u} \end{bmatrix} \tag{2.3}$$

where *Ut* represents the tangent input matrix of Equation (7).

The linearized system and input matrices can be computed for any working point at the time step *k* + 1 as:

$$A\_{k+1} = -G\_{\chi\_{k+1}}^{-1} G\_{\chi\_k};\tag{24}$$

$$B\_{k+1} = -G\_{x\_{k+1}}^{-1} G\_{\mu\_{k+1}}.\tag{25}$$

These equations enable the evaluation of the linearized time-discretized explicit description of the EOMs suitable for estimation methods such as the extended Kalman filter. Their deployment of this estimation scheme is discussed in the next section.

It is important to note that the invertibility of the matrix *Gxk*+1 is guaranteed thanks to the choice of BDF-1 differentiation scheme and that in contrast, it would not be possible. For completeness, the limitations of a forward differentiation scheme is demonstrated in the Appendix A.

#### **4. State-Input Estimation for MB Models**

The augmented discrete extended Kalman filter (ADE-KF) tailored for MB models is discussed in this section. In this section we discuss all the required components to set up a Kalman filter for assimilating the different multibody states and unknown inputs. More in details, in Section 4.1, the general form of the discrete-time system and measurement model equations are summarized. In Section 4.2, the measurement equations are augmented with the constraint equations *φ*, to deal with the constrained nature of the MB EOMs, leading to a constrained estimation problem. Moreover, in Section 4.3 the adopted approach to compute the linearized measurement matrices *C* and *D* are presented since they are required in the estimation framework and not directly available. Finally, in Section 4.4 the ADE-KF is assembled to deal with the estimation of states and unknown inputs, and the Kalman filter steps tailored for MB models are reported.

#### *4.1. Model and Measurement Equations with Uncertainty*

The system of EOMs in Equation (4) described in the previous section are assumed as deterministic. However, in practice, various noise sources will disturb the behavior of the system. For the dynamic model equations, the process noise vector *ν*∗ *x*,*k*+1 influences the response:

$$g\_d(\mathbf{x}\_{k+1'}, \mathbf{x}\_{k'} \mathbf{u}\_{k+1}) = \vec{v}\_{\mathbf{x}, k+1'} \tag{26}$$

where *ν*˜*x*,*k*+1 is associated with the noise term *νx*,*k*+1 described by the following equation:

$$\mathbf{x}\_{k+1} = \mathbf{x}\_{k+1}^\* + \boldsymbol{\nu}\_{\mathbf{x},k+1}.\tag{27}$$

Here, *x*∗ *k*+1 is the deterministically varying state vector while *νx*,*k*+1 is a zero-mean noise term with a (possibly time-varying) covariance *Qk*.

The estimation framework relies on combining the prior model information with measurements *y* ∈ R*ny* on the physical process. To compare the prediction of the model with these measurements, a set of measurement equations *h* are necessary which are affected by measurement noise *νy*,*k*+1:

$$y = h(\dot{\upsilon}\_{k+1}, \chi\_{k+1}, \mu\_{k+1}) + \upsilon\_{\circ, k+1}.\tag{28}$$

where again *νy*,*k*+1 is assumed to be zero-mean with a (time-varying) covariance *Rk*. Similar to the model measurement equations *y*, the concept of virtual sensor (VS) *yvs* ∈ R*nyvs* is introduced:

$$y\_{\text{vs}} = h\_{\text{vs}}(\dot{v}\_{k+1}, \mathbf{x}\_{k+1}, u\_{k+1}) + \boldsymbol{\upsilon}\_{\text{vs}}.\tag{29}$$

A virtual sensor represents model-based information; more specifically, it is a physical quantity that can be expressed as a function of the system variables, such as joint forces, body motion trajectories (e.g., linear and angular positions, velocities, and accelerations) and body strains/stresses. A VS equation can be treated similarly to a measurement equation but evaluated after the a posteriori Kalman filter step and hence to exploited in the estimation process itself.

#### *4.2. The Augmented Constraint Measurement Equations*

In [23] different approaches have been proposed to deal with constrained state estimation, although the authors deem it impossible to give a general guideline for constrained filtering because each individual problem is unique and dependent on how the constraints are handled.

In this work, the combination of hard and soft constraint methods is considered [23], where the constraint measurement equations are implemented in the Kalman filter scheme by augmenting the measurement equations *y* of Equation (28) with the constraint equations:

$$\Psi = \begin{cases} h(\dot{\nu}\_{k+1}, \mathbf{x}\_{k+1}, \mathbf{u}\_{k+1}) + \nu\_{\mathcal{Y}, k+1} \\ \phi(q\_{k+1}) + \nu\_{\phi, k+1} \end{cases} . \tag{30}$$

Here, *νφ*,*k*+1 is a small zero-mean noise term added to the idealized constraint equations *φ*(*qk*+1) if small constraints violation is allowed (soft constraints method). In the case *νφ*,*k*+1 is assumed to be zero, the constraint equations are treated as perfect measurements (hard constraint method).

#### *4.3. An Efficient Strategy for the Measurement Sensitivities Computation*

To evaluate the EKF equations at each filter step, the linearized measurement sensitivity matrices *Ck*+1 and *Dk*+1 are required and they are defined as:

$$\mathbf{C}\_{k+1} = \frac{\partial \overline{y}}{\partial \mathbf{x}}\Big|\_{k+1} = \begin{bmatrix} \frac{dy}{d\overline{x}}\\ \frac{d\Phi}{d\mathbf{x}} \end{bmatrix}\_{k+1} = \begin{bmatrix} \frac{dy}{d\overline{v}} & \frac{dy}{d\overline{q}} & \frac{dy}{d\overline{\Lambda}}\\ 0\_{\lambda,v} & J & 0 \end{bmatrix}\_{k+1};\tag{31}$$

$$D\_{k+1} = \frac{\partial \bar{y}}{\partial u}\Big|\_{k+1} = \begin{bmatrix} \frac{dy}{du} \\ 0\_{\lambda,u} \end{bmatrix}\_{k+1}.\tag{32}$$

Generally, the measurement equations can be expressed as a function of the generalized accelerations, velocities, positions, Lagrange multipliers, and inputs, Equation (30).

However, due to the state-space description presented in Section 3, the explicit relationship of the sensor equations with respect to the generalized accelerations is lost. More in particular, when considering acceleration sensors, their sensitivity with respect to the generalized accelerations is non zero, thus it must be indirectly included. Acceleration sensors are really common in the mechanical engineering community due to their noninvasive and simple installation requirements. Within the estimation framework, they can be used either in form of measurements or in form of VSs.

Similarly, the sensitivity of the measurement equations with respect to the system inputs are not always available or nontrivial to obtain.

Therefore, an approach to explicitly obtain these dependencies is here derived by means of the acceleration level constraints.

Introducing the auxiliary variables *z* ∈ R3*nq*+*nλ*

$$z = \begin{bmatrix} \dot{q}^T & \dot{q}^T & q^T & \lambda^T \end{bmatrix}^T = \begin{bmatrix} \dot{v}^T & v^T & q^T & \lambda^T \end{bmatrix}^T,\tag{33}$$

for a generic measurement sensor *y*, the sensitivities *∂y ∂z* with respect to the auxiliary variable *z* can be analytically computed starting from their fundamental equations.

$$\frac{\partial \underline{y}}{\partial z} = \begin{bmatrix} \frac{\partial y}{\partial v} & \frac{\partial y}{\partial v} & \frac{\partial y}{\partial q} & \frac{\partial y}{\partial \lambda} \end{bmatrix}. \tag{3.4}$$

Applying the chain rule, the non-augmented linear measurement matrix *dy dx* for a generic time step can be written as

$$\frac{dy}{dx} = \frac{\partial y}{\partial z}\frac{dz}{dx}.\tag{35}$$

And the full derivative *dz dx* can be expressed as

$$\frac{dz}{dx} = \begin{bmatrix} \frac{\partial \bar{x}}{\partial v} & \frac{\partial \bar{x}}{\partial q} & \frac{\partial \bar{x}}{\partial \lambda} \\ \frac{\partial v}{\partial v} & \frac{\partial v}{\partial q} & \frac{\partial v}{\partial \lambda} \\ \frac{\partial q}{\partial v} & \frac{\partial q}{\partial q} & \frac{\partial q}{\partial \lambda} \end{bmatrix} = \begin{bmatrix} \frac{\partial \psi}{\partial v} & \frac{\partial \psi}{\partial q} & \frac{\partial \psi}{\partial \lambda} \\ I\_v & 0\_{v,q} & 0\_{v,\lambda} \\ 0\_{q,v} & I\_q & 0\_{q,\lambda} \\ \frac{\partial \lambda}{\partial v} & \frac{\partial \bar{\lambda}}{\partial q} & I\_\lambda \end{bmatrix}. \tag{36}$$

To analytically compute the unknown terms of Equation (36), we introduce the acceleration level constraints *φ*¨ defined as the second time derivative of the position level constraints. Starting from the velocity level constraints *φ*˙:

$$
\dot{\phi} = \frac{\partial \phi}{\partial q} \frac{dq}{dt} = Jv,\tag{37}
$$

the acceleration level constraint equations are obtained as follows:

$$\ddot{\boldsymbol{\phi}} = \boldsymbol{J}\boldsymbol{\upsilon} + \boldsymbol{J}\dot{\boldsymbol{\upsilon}} = \sum\_{i,j=1}^{n\_q} \boldsymbol{\upsilon}\_i \frac{\partial^2 \boldsymbol{\phi}}{\partial [q\_i, q\_j]} \boldsymbol{\upsilon}\_j + \boldsymbol{J}\dot{\boldsymbol{\upsilon}} = \boldsymbol{0}\_\lambda. \tag{38}$$

By means of the continuous dynamic equations *g*2 of Equation (8) and the acceleration level constraint of Equation (38), a new set of equations *p* ∈ R*nq*+*nλ* can be constructed:

$$p(z,u) = \begin{cases} p\_1 = M(q)\dot{\upsilon} + f\_{nl}(\upsilon, q\_\prime u) + f^T \lambda = 0\_\emptyset \\ p\_2 = \sum\_{i,j=1}^{n\_\emptyset} \upsilon\_i \frac{\partial^2 \phi}{\partial [q\_{i\prime} q\_j]} \upsilon\_j + f \dot{\upsilon} = 0\_\lambda \end{cases} \tag{39}$$

Equations (39) represent the governing set of equations that implicitly determines the relation between the different coordinates *v*, *q*, *v*˙ and the Lagrange multipliers *λ*.

Therefore, the unknown terms of Equation (36) can be computed as

$$
\begin{split}
\begin{bmatrix}
\frac{\partial\psi}{\partial\upsilon} & \frac{\partial\psi}{\partial q} & \frac{\partial\psi}{\partial\lambda} \\
\frac{\partial\lambda}{\partial\upsilon} & \frac{\partial\lambda}{\partial q} & \frac{\partial\lambda}{\partial\lambda}
\end{bmatrix} &= -\begin{bmatrix}
\frac{\partial p}{\partial\upsilon} & \frac{\partial p}{\partial\lambda}
\end{bmatrix}^{-1} \begin{bmatrix}
\frac{\partial p}{\partial\upsilon} & \frac{\partial p}{\partial q} & \frac{\partial p}{\partial\lambda}
\end{bmatrix} = \\
& -\begin{bmatrix}
\mathcal{M}(q) & f^{T} \\
f & 0\_{\lambda,\lambda}
\end{bmatrix}^{-1} \begin{bmatrix}
\mathcal{C}\_{t} & \mathcal{K}\_{t} & f^{T} \\
2\sum\_{j=1}^{n\_{\upsilon}} \frac{\partial^{2}\phi}{\partial[q\_{j}q\_{j}]} \upsilon\_{j} & \sum\_{j=1}^{n\_{\upsilon}} \frac{\partial^{2}\phi}{\partial[q\_{j}q\_{j}]} \upsilon\_{j} & 0\_{\lambda,\lambda}
\end{bmatrix}
\end{split}
\tag{40}
$$

where the third order partial derivative of the constraint equations *φ*(*q*) with respect to the generalized coordinates *q* is zero for FNCF. Another important ingredient to fulfill the Kalman filter requirements is the measurement input matrix *dy du* of Equation (32). However, this matrix is not directly available and it is computed similarly to Equation (35) but with respect to input *u*:

$$\frac{d\underline{y}}{d\underline{u}} = \frac{\partial y}{\partial z}\frac{dz}{d\underline{u}} + \frac{\partial \underline{y}}{\partial \underline{u}}.\tag{41}$$

Here, the full derivative *dz du* represents how the vector *z* varies when the input *u* is perturbed obtaining:

$$\frac{dz}{du} = \begin{bmatrix} \frac{\partial \psi}{\partial u} & \frac{\partial \psi}{\partial u} & \frac{\partial q}{\partial u} & \frac{\partial \lambda}{\partial u} \end{bmatrix} = \begin{bmatrix} \frac{\partial \psi}{\partial u} & 0\_{\upsilon,u} & 0\_{q,u} & \frac{\partial \lambda}{\partial u} \end{bmatrix}.\tag{42}$$

From Equation (42) can be seen that only the sensitivities of the acceleration and Lagrange multipliers are non-zero terms. A second order system is fully defined by the position coordinates, velocity coordinates, and time. An external force only influences the force balance of that system and thus the acceleration coordinates and Lagrangian multipliers while only through time-integration the velocity and position coordinates. However, these derivatives are evaluated at a single time instance. Therefore, in Equation (36), the acceleration coordinates and Lagrangian multipliers depend on the position and velocity coordinates, but not the other way around: the position and velocity coordinates do not depend on the acceleration coordinates, Lagrangian multipliers, or the external forces.

Subsequently, the remaining terms of Equation (42) can be computed as

$$
\begin{bmatrix}
\frac{\partial \psi}{\partial \lambda} \\
\frac{\partial \lambda}{\partial \mu}
\end{bmatrix} = -\begin{bmatrix}
\frac{\partial p}{\partial \psi} & \frac{\partial p}{\partial \lambda}
\end{bmatrix}^{-1} \frac{\partial p}{\partial \mu} = -\begin{bmatrix}
M(q) & J^T \\
J & 0\_{\lambda,\lambda}
\end{bmatrix}^{-1} \begin{bmatrix}
\mathcal{U}\_{\mathbb{I}} \\
0\_{\lambda,\mu}
\end{bmatrix}.\tag{43}
$$

In Figure 1 a schematic summary of the followed steps required for the computation of the system and measurement matrices is reported.

**Figure 1.** Block diagram representation of the system and measurement matrices computation for a generic integration step.

#### *4.4. Augmented Discrete Extended Kalman Filter*

The goal of this work is to estimate the unknown input forces, often referred to as disturbances, *d* ∈ R*nd* , concurrently with the system states *x* through the augmented state vector *x*˜ ∈ R*nx*+*nd* defined as:

$$\mathfrak{x}\_{k+1} = \begin{bmatrix} \mathfrak{x}\_{k+1} \\ d\_{k+1} \end{bmatrix}. \tag{44}$$

As the estimation algorithm also needs a model for the unknown inputs, a random walk model is associated with the unknown input:

$$d\_{k+1} = d\_k + \nu\_{d,k};\tag{45}$$

with *νd*,*k* also being zero-mean, random noise with covariance matrix *Qd*,*k*.

In this work, the random walk model is chosen for the unknown disturbance. This approach is chosen as in general no prior information on the input is assumed to be known. However, if useful information is available, e.g., periodicity of the input, other input models can be employed, as in [24].

The augmented discrete-time state-space system of equations can be written combining Equations (26) and (45):

$$\mathfrak{g}\_d(\mathfrak{x}\_{k+1}, \mathfrak{x}\_k, \mathfrak{u}\_{k+1}) = \begin{cases} \mathfrak{g}\_d(\mathfrak{x}\_{k+1}, \mathfrak{x}\_{k'} \mathfrak{u}\_{k+1}) = \mathfrak{v}\_{\mathfrak{x},k} \\ d\_{k+1} - d\_k = \nu\_{d,k} \end{cases} . \tag{46}$$

The linearized system matrix F refers to the augmented system of equations Equation (46) at the time step *k* + 1 and it is assembled as

$$\mathcal{F} = \frac{\partial \overline{f}\_d}{\partial \overline{\mathfrak{X}}} \Big|\_{k+1} = \begin{bmatrix} A\_{k+1} & B\_{k+1}^d \\ 0\_{\mathbf{x},d} & I\_d \end{bmatrix};\tag{47}$$

where ˜ *fd* is the augmented explicit ODEs associated with Equation (46), is obtained from Equation (24), while *Bd k*+1 is the linear disturbance matrix obtained similarly to *Bk*+1 in Equation (25) but with respect to the unknown input *d*.

Subsequently, the linearized augmented measurement matrix H associated with the augmented state vector *x*˜ is assembled as:

$$\mathcal{H} = \frac{\partial \vec{y}}{\partial \vec{x}}\Big|\_{k+1} = \begin{bmatrix} \mathbb{C}\_{k+1} & D\_{k+1}^d \end{bmatrix} \tag{48}$$

where the matrix *Dd k*+1 is computed in a similar fashion as for *Dk*+1 in Equation (32) but with respect to the unknown input *d*.

Starting from these model, measurement and tangent matrices, an extended Kalman filter (EKF) procedure can be set up in order to estimate both the states and the unknown input forces. The application of this estimation approach is discussed in the following section.

#### *4.5. The Adopted Extended Kalman Filter Scheme*

In this work, we employ an extended Kalman filter (EKF) [25] to perform the stateestimation. With the various model and measurement equations, and their respective derivatives, defined for a MB model in the previous sections this EKF can be run efficiently. As all terms have been defined analytically it is not necessary to resort to numerical derivatives, which allows for an effective application of the EKF. Several researchers have shown the applicability of the estimators even for the strongly nonlinear case of MB systems [9,12,26].

In general, we can group the EKF scheme in two main steps: the a-priori and a posteriori steps.

• *A-priori step*: Assuming that the augmented states *x*˜ + *k* at the previous filter step and the input *uk*+1 are known, the a-priori state prediction *x*˜ − *k*+1 and generalized accelerations *v*˙ − *k*+1 can be computed solving the ID-DAEs of Equation (17):

$$\mathbf{g}\_d(\tilde{\mathbf{x}}\_{k+1'}^-, \tilde{\mathbf{x}}\_{k'}^+, \boldsymbol{\mu}\_{k+1}) = \mathbf{0}\_{\mathbf{x}}.\tag{49}$$

Knowing the estimated state covariance matrix P+ *k* for the previous timestep, the apriori covariance at the current time (*k* +1) step can be approximated from Equation (47) as

$$\mathcal{P}\_{k+1}^{-} = \mathcal{F}\mathcal{P}\_{k}^{+}\mathcal{F}^{T} + \bar{\mathcal{Q}}\_{k}.\tag{50}$$

where *Q*˜ *k* is the process noise matrix of the augmented state-space model. The predicted measurement *y*˜ − *k*+1 can then be evaluated from Equation (30) as:

$$
\mathfrak{Y}\_{k+1}^- = \mathfrak{Y}(\mathfrak{v}\_{k+1'}^-, \mathfrak{x}\_{k+1'}^-, \mathfrak{u}\_{k+1}).\tag{51}
$$

The Kalman filter gain K allows achieving a desireable trade-off between the confidence in the model and the available measurements, and can be evaluated as:

$$\mathcal{K}\_{k+1} = \mathcal{P}\_{k+1}^{-} \mathcal{H}^{T} (\mathcal{H} \mathcal{P}\_{k+1}^{-} \mathcal{H}^{T} + \mathbb{R}\_{k})^{-1},\tag{52}$$

where H is obtained from Equation (48) and *R*˜ *k* is the measurement covariance matrix of the augmented measurement equations.

• *A-posteriori step*: When the real measurement *y*∗ *k*+1 becomes available together with the predicted measurement *y*˜ − *k*+1, the a posteriori state vector *x*˜ + *k*+1 is obtained as:

$$
\mathfrak{x}\_{k+1}^{+} = \mathfrak{x}\_{k+1}^{-} + \mathcal{K}\_{k+1} (y\_{k+1}^{\*} - \mathfrak{y}\_{k+1}^{-}) \tag{53}
$$

The inclusion of the actual measurements also affects the posterior covariance matrix P+ *k*+1 and can be evaluated as:

$$
\mathcal{P}\_{k+1}^{+} = (\mathcal{I}\_{\mathbf{x}} - \mathcal{K}\_{k+1} \mathcal{H}) \mathcal{P}\_{k+1}^{-}. \tag{54}
$$

In Figure 2 a block diagram scheme summarizing a generic *k*th step of the recursive ADE-KF is shown.

**Figure 2.** Block diagram representation of the recursive ADE-KF scheme for a generic filter step.

With this scheme the entire state-input estimation can be performed for an arbitrary (flexible) multibody model. In the following section we will validate it on an experimental setup.

In practical applications the performance of the estimation scheme will strongly depend on the tuning of the model covariance *Q*˜ and measurement covariance *R*˜. Unfortunately this tuning is also highly case dependent, which makes it difficult to set up general guidelines. For the particular validation case considered in this work, a detailed discussion on a possible tuning strategy is presented in the following section.

#### **5. Validation: Joint State-Input Estimation**

In this section, the presented joint state-input estimation methodology is experimentally validated on a slider-crank mechanism.

#### *5.1. The Slider-Crank System*

The slider-crank system in Figure 3 is a mechanism widely used in industry to transform a rotational motion in a linear motion as for instance in pick and place operations or in car engines.

**Figure 3.** The slider-crank: experimental setup.

Figure 4 represents the multibody model of the experimental setup shown in Figure 3 consisting of 14 bodies: base block (BB), motor housing (*MH*), motor rotor (*MR*), motor support (*MS*), crank (*C*), crank shaft (*CSh*), crank support (*CS*), connecting rod (*CR*), crankrod bearing (*BC*−*CR*), slider (*S*), rod-slider bearing (*BCR*−*S*), slider accelerometer (*SA*), track rail (*TR*) and track support (*TS*).

**Figure 4.** The slider-crank: multibody model.

The mechanical properties of the various bodies are listed in Table 1.

**Table 1.** Mechanical properties expressed with respect to the center of gravity of each individual body.

Even though the presented slider-crank mechanism is an academic example it combines several challenging phenomena such as non-linear kinematics and complex slidertrack interaction phenomena. The various bodies are connected trough fixed, revolute, spherical, and universal joints to allow the desired kinematics. The slider and track bodies are connected through a contact stiffness **k***c* along the perpendicular directions to the sliding trajectory. The Pacejka friction model [27] defines the friction coefficient *μ* as a function of the sliding velocity Δ*v*:

$$\mu = \mathbf{d} \sin \left[ \mathfrak{c} \tan^{-1} \left( \mathbf{b} \Delta v - \mathbf{e} \left( \mathbf{b} \Delta v - \tan^{-1} (\mathbf{b} \Delta v) \right) \right) \right],\tag{55}$$

where the friction force *f f* is evaluated as

$$f\_f = \mu f\_\mathbb{n} = \mu(\mathbf{k}\_\mathbb{c}\delta + \mathbf{c}\_\mathbb{c}\dot{\delta}),\tag{56}$$

where *fn* is the resultant normal force acting on the slider and, *δ* is the local compliance between the slider and track rail interfaces projected onto the normal direction to the sliding direction. The adopted (identified) Pacejka friction model parameters are listed in Table 2.

**Table 2.** Pacejka friction model parameters.

The system motion is driven by a brushless servomotor "MAC3000"' with integrated controller MAC00-B4 from JVL (www.jvl.dk, accessed on 1 March 2021). It includes a motor encoder sensor allowing the measurement of the rotation angle and velocity of the motor shaft, indicated by *θ* and ˙ *θ* respectively.As can be seen in Figures 3 and 4, the slider is equipped with a mono-axial MEMS accelerometer (3711D1FB200G) from PCB (www.pcb.com) to measure its translational acceleration along the sliding direction indicated by *Y*¨. The entire system is controlled via the motor using a closed-loop PID controller targeting a desired rotational motor velocity while adapting the determined motor torque *T*.


**Meet the editor**

Anna Rudawska, PhD, DSc, Eng, Assoc Prof, is the head of the Technological Processes and Systems Design Chair, and since September 1, 2012, she has been the vice dean of the Faculty of Mechanical Engineering at Lublin University of Technology, Poland. She has authored or coauthored 200 scientific publications and several patents in the adhesive field. She is a member of many

international and national organizations and a reviewer and a member of editorial boards of international journals. Her research focuses on the issues of analysis of the influence of technological and structural factors on the adhesive joint strength of polymers and metals, design of bonding technology, adhesive properties research, as well as issues related to receiving the appropriate adhesive properties to increase the strength of the adhesive joint.

Contents

**Preface VII**

Jingsi Wang

**Modeling 65**

**Abrasive Methods 1**

**Inside AWJ Nozzle 3** Xinping Long and Qi Liu

**Section 1 Characteristics and Application of Abrasive Technology and**

Chapter 1 **Numerical Analysis of the Liquid-Gas-Solid Three Phase Flow**

Chapter 3 **Advanced Abrasive Waterjet for Multimode Machining 39**

Chapter 5 **HyDRa: Vortex Polishing with a Deterministic Hydrodynamic**

**Section 2 Characteristics of Abrasive Tools, Materials and Machines 105**

H.-T. (Peter) Liu, Vanessa Cutler, Chidambaram Raghavan, Peter

Erika Sohn, Esteban Luna, Elfego Ruiz, Luis Salas and Joel Herrera

Artur Jamrozik, Łukasz Klapiszewski, Beata Strzemiecka, Adam

Chapter 2 **Ultrasonic Machining: A Total Mechanical Machining Technology Using Loose Abrasive Particles 19**

Miles, Ernst Schubert and Nathan Webers

Pavel Kovač and Marin Gostimirović

**Radial Polishing Tool 83**

Chapter 6 **Additives for Abrasive Materials 107**

Voelkel and Teofil Jesionowski

Chapter 4 **Grinding Force of Cylindrical and Creep-Feed Grinding**

## Contents

**Preface XI**

Preface

*als, and Machines.*

chining, and laser micro-machining technologies.

and their characteristics, and the type of machining materials.

The aim of this book is to present information on the characteristics and applications of abra‐ sive technology, abrasive tools, and tests, as well as the innovative methods of this technolo‐ gy. This book presents information regarding the abrasive technology methods, the characteristics of the methods (for example, the technological parameters, tools, and ma‐ chines), innovative methods, characteristics of surface structure and surface properties after this type of mechanical process, the tests of various issues of abrasive technology, and also the application in various industrial branches and other technical and technological domains. Abrasive technology is becoming increasingly important in precision manufacturing. It cov‐ ers the mechanics and mechanisms of abrasive processes, characteristics of abrasive tools, their structure, materials, and wear, as well as the grinding, finishing, sanding, deburring, hybrid super-finishing, and the abrasive-jet and flow machining processes, among others. Moreover, this book includes the issues on high-speed machining, eco-machining, nano-ma‐

There are a lot of factors, which influence the abrasive process—abrasive tools, technological parameters of the abrasive process, the type of machining methods, the type of machines

One of the important factors, which has a significant influence on the abrasive process and the surface roughness, is the type, structure, and properties of the abrasive material. Abra‐ sives generally rely on a difference in hardness between the abrasive and the material being worked upon, the abrasive being the harder of the two substances. The shape, size, and na‐ ture of the workpiece and the desired finish will influence the choice of the abrasive used. For example, a much harder abrasive will cut faster and deeper. The adhesion between grains, between grains and backing, between grains and matrix determines how quickly grains are lost from the abrasive and how soon fresh grains, if present, are exposed. More‐ over, there are lots of abrasive materials. In the past, abrasive grains were derived from nat‐ ural minerals, such as emery, quartz sand, and so on. Today, the industry is almost entirely

dependent upon synthetic grains, which are far more consistent and cost-effective.

branches of industry, as well as in dentistry and medical applications.

The book was divided into two parts: the first part: *Characteristics and Applications of Abrasive Technology and Abrasive Methods* and the second part: *Characteristics of Abrasive Tools, Materi‐*

The aim of such a presentation was to present the usage of characteristics and application of the abrasive methods (also innovative methods and tests) as well as the abrasive tools in very different aspects. These issues are very important in precision manufacturing in many

Chapter 10 **Abrasive for Chemical Mechanical Polishing 183** Hong Jin Kim

## Preface

Chapter 7 **Tuning the Morphology and Surface Property of Mineral**

Chapter 8 **Recycling and Recharging of Supreme Garnet in Abrasive**

Chapter 9 **Nanostructural Deformation in Brittle-Ductile Compounds and Its Application in Fabricating Ductile Nanoparticles 161**

**Particles by Grinding Media 127** Zhiyong Gao and Chengwei Li

**Waterjet Machining 147** Vu Ngoc Pi and Tran Quoc Hung

Babak Alinejad and Yazdan Zare

Hong Jin Kim

**VI** Contents

Chapter 10 **Abrasive for Chemical Mechanical Polishing 183**

The aim of this book is to present information on the characteristics and applications of abra‐ sive technology, abrasive tools, and tests, as well as the innovative methods of this technolo‐ gy. This book presents information regarding the abrasive technology methods, the characteristics of the methods (for example, the technological parameters, tools, and ma‐ chines), innovative methods, characteristics of surface structure and surface properties after this type of mechanical process, the tests of various issues of abrasive technology, and also the application in various industrial branches and other technical and technological domains.

Abrasive technology is becoming increasingly important in precision manufacturing. It cov‐ ers the mechanics and mechanisms of abrasive processes, characteristics of abrasive tools, their structure, materials, and wear, as well as the grinding, finishing, sanding, deburring, hybrid super-finishing, and the abrasive-jet and flow machining processes, among others. Moreover, this book includes the issues on high-speed machining, eco-machining, nano-ma‐ chining, and laser micro-machining technologies.

There are a lot of factors, which influence the abrasive process—abrasive tools, technological parameters of the abrasive process, the type of machining methods, the type of machines and their characteristics, and the type of machining materials.

One of the important factors, which has a significant influence on the abrasive process and the surface roughness, is the type, structure, and properties of the abrasive material. Abra‐ sives generally rely on a difference in hardness between the abrasive and the material being worked upon, the abrasive being the harder of the two substances. The shape, size, and na‐ ture of the workpiece and the desired finish will influence the choice of the abrasive used. For example, a much harder abrasive will cut faster and deeper. The adhesion between grains, between grains and backing, between grains and matrix determines how quickly grains are lost from the abrasive and how soon fresh grains, if present, are exposed. More‐ over, there are lots of abrasive materials. In the past, abrasive grains were derived from nat‐ ural minerals, such as emery, quartz sand, and so on. Today, the industry is almost entirely dependent upon synthetic grains, which are far more consistent and cost-effective.

The book was divided into two parts: the first part: *Characteristics and Applications of Abrasive Technology and Abrasive Methods* and the second part: *Characteristics of Abrasive Tools, Materi‐ als, and Machines.*

The aim of such a presentation was to present the usage of characteristics and application of the abrasive methods (also innovative methods and tests) as well as the abrasive tools in very different aspects. These issues are very important in precision manufacturing in many branches of industry, as well as in dentistry and medical applications.

This information enables scientists, engineers, and designers to ensure the soundness and integrity of the fabricated components and to develop new techniques effectively. I would like to express my gratitude to all the contributors for high-quality manuscripts. I hope that the open-access format of this book will help researchers and engineers.

**Anna Rudawska**

**Section 1**

**Characteristics and Application of Abrasive**

**Technology and Abrasive Methods**

Lublin University of Technology Faculty of Mechanical Engineering Department of Production Engineering Poland