Robert J Riggs
email: robriggs@umich.edu

Olga Battaïa

S Jack Hu

Disassembly line balancing under high variety of end of life states using a joint precedence graph approach

Disassembly line balancing (DLB) Joint precedence graph End of life (EOL) condition Stochastic line balancing a b s t r a c t

Disassembly is an important aspect of end of life product treatment, as well as having products disassembled in an efficient and responsible manner. Disassembly line balancing is a technique that enables a product to be disassembled as efficiently and economically viable as possible; however, considering all possible end of life (EOL) states of a product makes disassembly line balancing very difficult. The EOL state and the possibility of multiple recovery options of a product can alter both disassembly tasks and task times for the disassembly of the EOL product. This paper shows how generating a joint precedence graph based on the different EOL states of a product is beneficial to achieving an optimal line balance where traditional line balancing approaches are used. We use a simple example of a pen from the literature to show how a joint disassembly precedence graph is created and a laptop example for joint precedence graph generation and balancing. We run multiple scenarios where the EOL conditions have different probabilities and compare results for the case of deterministic task times. We also consider the possibility where some disassembly task times are normally distributed and show how a stochastic joint precedence graph can be created and used in a stochastic line balancing formulation.

Introduction

The end of life treatment of products is a very important topic concerning manufacturers, consumers, governments, and society as a whole. The waste from these products can have negative impacts on the environment creating a less sustainable future. Many national governments (e.g., Japan, Canada, and Taiwan), the European Union members, and 23 states in the US have adopted extended producer responsibility (EPR) principle based legislation for end-of-use treatment of products [START_REF] Ozdemir | Recovery decisions of a producer in a legislative disposal fee environment[END_REF]. Due to this fact, manufacturers are increasingly recovering, remanufacturing, and reprocessing postconsumer products. For instance, the European Union (EU) created directives requiring companies to be responsible, free of charge, for the end of life treatment of many products including electronics and automobiles. The Waste Electrical and Electronic Equipment (WEEE) directive, together with the Restriction of Hazardous Substances (RoHS) directive, imposes the responsibility of disposal of electrical waste and electronic equipment to the manufacturers of the equipment. Disassembly in some form is the common thread for all forms of product recovery. Product recovery required by law is being discussed and implemented around the world as countries become conscious of the need for sustainable practices and the potential economic benefit of remanufacturing.

Developing a disassembly system enables products to be disassembled in a responsible and efficient manner. There are many decisions to be made for the design and layout of a disassembly system, such as the number of workstations and task assignment. Line balancing (LB) is a decision making tool that assigns manufacturing tasks to a set number of workstations with the objective of minimizing some performance objective, such as the maximum difference in total task time between workstations. LB was first considered in assembly settings but the approach can be extended to the application of disassembly [START_REF] Mcgovern | The disassembly line: balancing and modeling[END_REF][START_REF] Battaïa | A taxonomy of line balancing problems and their solution approaches[END_REF]. Meacham et al. [START_REF] Meacham | Optimal disassembly configurations for single and multiple products[END_REF] determined the optimal disassembly configurations for both single and multiple product types for meeting a specified demand for recovered components and subassemblies. Gungor and Gupta [START_REF] Gungor | A solution approach to the disassembly line balancing problem in the presence of task failures[END_REF] investigated the disassembly line balancing problem (DLBP) in the presence of task failures. Some disassembly tasks may not be completed due to a product defect that affects the rest of the disassembly process and some disassembly steps may need to be altered or skipped entirely. Gungor et al. [START_REF] Gungor | Complications in disassembly line balancing[END_REF] also discussed the challenges that come from disassembly line balancing, ranging from product, disassembly line, demand, and assignment complications.

Altekin et al. [START_REF] Altekin | Disassembly line balancing with limited supply and subassembly availability[END_REF] investigated the DLBP as a partial disassembly with limited supply of a single product and its subassemblies. They created two formulations, one that maximizes profit over a single disassembly cycle and the other maximizes profit over the entire planning horizon. Tang and Zhou [START_REF] Tang | A systematic approach to design and operation of disassembly lines[END_REF] created a disassembly Petri net model for the hierarchical modeling in order to derive a disassembly path with the maximal benefit in the presence of some defective components. Their algorithm for balancing disassembly lines seeks to maximize the productivity of a disassembly system. McGovern et al. [START_REF] Mcgovern | A balancing method and genetic algorithm for disassembly line balancing[END_REF] proved that the disassembly line balancing problem is NP-complete and presents a genetic algorithm for balancing disassembly lines and used priori instances to evaluate any disassembly line balancing technique.

Altekin et al. [START_REF] Altekin | Profit-oriented disassembly-line balancing[END_REF] developed an MIP formulation for solving partial disassembly-line balancing that determines simultaneously multiple decision points such as the cycle time and task assignment of the line and the number of stations to be opened. Tripathi et al. [START_REF] Tripathi | Real world disassembly modeling and sequencing problem: Optimization by Algorithm of Self-Guided Ants (ASGA)[END_REF] proposed a fuzzy disassembly optimization model that determines the disassembly sequence and the depth of disassembly to maximize net revenue from EOL disposal of products. Fuzzy control theory is used to account for the uncertainty associated with the quality in returns. Koc et al. [START_REF] Koc | Two exact formulations for disassembly line balancing problems with task precedence diagram construction using an AND/OR graph[END_REF] developed both an integer and dynamic programming formulations for the DLBP by using AND/OR graphs to check for feasibility of the solution.

Altekin et al. [START_REF] Altekin | Task-failure-driven rebalancing of disassembly lines[END_REF] developed a two-step approach to rebalance a disassembly line when task failures occur that lead to successor tasks being infeasible. Ma et al. [START_REF] Ma | Disassembly process planning algorithms for end-of-life product recovery and environmentally conscious disposal[END_REF] created a model that simultaneously solves the decisions of disassembly level, sequence, and EOL options for components or subassemblies in a product in the parallel disassembly environment for the objective of maximizing profit. The parallel disassembly environment allows two or more parts of subassemblies to be disassembled at the same time. McGovern et al. [START_REF] Mcgovern | The disassembly line: balancing and modeling[END_REF] formulated the disassembly line balancing problem and presented case studies and solution methods in their book The Disassembly Line: Balancing and Modeling. Rickli et al. [START_REF] Rickli | Multi-objective partial disassembly optimization based on sequence feasibility[END_REF] created a multi-objective genetic algorithm that chose the optimal partial disassembly sequence with respect to operation and recovery costs, as well as revenue and environmental impact. Kalayci and Gupta [START_REF] Kalayci | Artificial bee colony algorithm for solving sequence-dependent disassembly line balancing problem[END_REF] used an artificial bee colony optimization technique to solve the disassembly line balancing problem with sequence dependence. They compare their bee colony algorithm to six metaheuristic approaches from the literature and demonstrate the effectiveness of the algorithm in solving the disassembly line balancing problem with multi-objectives. Riggs et al. [START_REF] Riggs | Two-stage sequence generation for partial disassembly of products with sequence dependent task times[END_REF] developed a two-stage sequence generation approach for the partial disassembly of products that have sequence dependent task times.

Ilgin and Gupta [START_REF] Ilgin | Environmentally conscious manufacturing and product recovery (ECMPRO): a review of the state of the art[END_REF] provided an extensive review of the state of the art for EOL product recovery and its relationship with disassembly for remanufacturing and recycling. Tuncel et al. [START_REF] Tuncel | Solving large scale disassembly line balancing problem with uncertainty using reinforcement learning[END_REF] used a reinforcement learning approach for disassembly line balancing where uncertainty comes from a single source, demand fluctuations for the disassembled components. They assumed all incoming products have the same identical structure and are in the same EOL condition. Other approaches have focused on uncertainty coming from the EOL product itself, which results in task failures or additional processing [START_REF] Gungor | A solution approach to the disassembly line balancing problem in the presence of task failures[END_REF][START_REF] Altekin | Task-failure-driven rebalancing of disassembly lines[END_REF].

There are many different approaches used in the literature to represent the relative order of disassembly of a product at its EOL. The approaches in the literature are referred to as either tree-based or state-based for the structure of disassembly [START_REF] Lambert | Disassembly sequencing: a survey[END_REF]. A tree-based approach assumes there is only one way to disassemble a product into its components and subassemblies, similar to a precedence graph. The state-based approach shows multiple possible disassembly paths, such as shown in an AND/OR graph or transition matrix [START_REF] Kang | Automatic derivation of transition matrix for end-of-life decision making[END_REF]. Many approaches in the literature use AND/OR graphs [START_REF] Aydemir-Karadag | Multi-objective optimization of stochastic disassembly line balancing with station paralleling[END_REF] to show the possibility of different disassembly paths; however, when using an AND/OR graph, only one possible disassembly path is chosen and the line balance resulting from that path must be used for all products at their EOL.

For stochastic disassembly line balancing, the literature has less breadth than the literature for deterministic disassembly line balancing. Agrawal et al. [START_REF] Agrawal | A collaborative ant colony algorithm to stochastic mixedmodel U-shaped disassembly line balancing and sequencing problem[END_REF] developed an ant colony algorithm to solve for a mixed-model U-shaped disassembly line to handle the situation where completion times are stochastic and there is task time variability due to the human factor. Aydemir-Karadag and Turkbey [START_REF] Aydemir-Karadag | Multi-objective optimization of stochastic disassembly line balancing with station paralleling[END_REF] developed a multi-objective genetic algorithm to determine the best line balance that considered stochastic task times where station paralleling is allowed. Bentaha et al. [START_REF] Bentaha | Dealing with uncertainty in disassembly line design[END_REF][START_REF] Bentaha | A sample average approximation method for disassembly line balancing problem under uncertainty[END_REF][START_REF] Bentaha | An exact solution approach for disassembly line balancing problem under uncertainty of the task processing times[END_REF] developed methods for disassembly line balancing where there is uncertainty in task processing times.

Many line balancing approaches for disassembly are the same as those used for assembly but with a simple inverse or alteration of the precedence graph [START_REF] Battaïa | A taxonomy of line balancing problems and their solution approaches[END_REF]. "Disassembly modeling and planning can be more challenging than assembly because its terminal goal is not necessarily fixed" due to changing market demand for reused components and subassemblies [START_REF] Tang | Disassembly modeling: planning, and application[END_REF]. The research for the disassembly line balancing problem (DLBP), up to this point, has not considered how the end of life (EOL) condition and treatment options for the same components will alter the precedence graph, either by representation or the time values for each task. At a product's EOL, disassembly of certain components or modules may have a different time distribution for the same disassembly task due to the variation in the EOL states. Some components can be grouped together and the entire module can be reused so complete disassembly is not required for each component. This will lead to tasks being skipped and having zero time. In addition, even if all EOL states are considered for the line balancing algorithm, it is difficult to include all this information and have the algorithm solve quickly. Previous approaches use AND/OR graphs to model the possibility of having different disassembly paths, but only one path is chosen and every product at its EOL takes this path for disassembly. Our approach is able to consider many disassembly paths that are weighted into a single graph.

Our objective in this paper is to develop and validate a method originally used in the area of mixed model assembly line balancing, joint precedence graph generation, and enhance it for the application of disassembly with multiple component/module EOL states. We treat each EOL state of a product as a different model variant and generate a disassembly joint precedence graph from all EOL states that can then be inputted into a line balancing algorithm. The contribution of our work is the development of a method that: (1) considers the existence of all EOL states with certain probabilities so that different line balancing techniques where a single precedence graph is required can be applied, (2) considers components/modules having multiple EOL treatment possibilities, and (3) is tractable and reduces the complexity of dealing with many EOL states and treatment possibilities during the optimization stage (i.e. we preprocess the complexity of many EOL states during the formulation of the joint precedence graph). To the best of our knowledge, this is the first method for balancing disassembly lines that is able to consider many EOL states and treatment options into a single precedence graph.

The remainder of the paper is organized as follows: Section 2 presents the method of generating disassembly joint precedence graphs and stochastic disassembly joint precedence graph. Section 3 outlines the deterministic and stochastic line balancing formulations. Section 4 presents line balancing results from a theoretical example for disassembly of a laptop taken from the literature and Section 5 contains discussion. Conclusions are discussed in Section 6.

Generation of disassembly joint precedence graphs

A precedence graph is a type of directed graph where arcs connect nodes to one another. The nodes are the tasks required for disassembly and the connecting arcs have a particular direction that shows the order of assembly/disassembly tasks. Boysen et al. [START_REF] Boysen | Assembly line balancing: joint precedence graphs under high product variety[END_REF] use a joint precedence graph for balancing mixed-model assembly lines by considering all products to be assembled on the same assembly line. This approach can be used to model precedence relations for disassembly lines where instead of considering different models of a product, the joint precedence graph will reflect every possible EOL state of a product. The disassembly line can then be balanced using existing line balancing optimization techniques from the literature. The EOL condition refers to the condition of an individual component or module and the EOL state is the combination of all EOL conditions for all components and modules in the assembly. The EOL state will have a probability associated with it and the summation of all EOL states will add up to one.

Method for generating joint precedence graphs

The starting point for creating a disassembly joint precedence graph for a product is to create a baseline disassembly precedence graph based on the product having "off the assembly line" quality. A product with off the assembly line quality has no quality defects and the EOL condition is like new. The baseline disassembly precedence graph should reflect all disassembly tasks that need to be accomplished to have the product disassembled to the required ending state.

Through various sampling methods, products at their end of life can be collected and the impact of different types of damage to the core can be assessed for how it can affect disassembly task times and the probability that type of damage will present itself. Damage to a product can increase or decrease task times resulting in variation of individual task times. Damage may also cause multiple tasks times to change simultaneously. For example, there is a disassembly task to remove component X and a component Y, first X and then Y. Component X has a 25% chance of being damaged, and the task time to remove component X is greater than the baseline case if component X is damaged. If component X is damaged, then the chance that component Y is damaged is 50% and the task time to remove Y is also greater than the baseline case, but if component X is not damaged then component Y would not be damaged and the baseline case for each task time is unchanged. This is an example of a dependency between two tasks. An example of independence between disassembly tasks is that damage to component X has no effect on the condition of component Y. Another possible phenomenon when disassembling a product is a task failure due to either damage to the product or absence of a component or module. Certain components or modules may become loose or dislodged during transport to its EOL treatment site or consumers may remove certain components or modules before disposing of the product. Both these cases can cause certain disassembly task times to be zero or negligible.

Once all possible component and module EOL states are determined, including both the difference in disassembly task times and the probability of having a difference in certain task times versus the baseline case, through simple combinatorics the probability of each EOL state can be enumerated and subsequent EOL state disassembly precedence graphs created. Since every state is enumerated, no two states will have exactly the same EOL profile. The EOL state disassembly precedence graphs will be very similar to the baseline disassembly precedence graph, except now a probability will be associated with each precedence graph, some tasks may be skipped, and the task times may be different for certain tasks.

The joint precedence graph is created by a weighted average approach that considers all precedence graphs for each EOL state. Task times in each precedence graph will be referred to as t iq , where i is the task number and q is the EOL state. Each EOL state will have a probability p q , where 0 ≤ p q ≤ 1, and the sum of all p q values for q Q is equal to 1, Q being the set of all EOL states. The average task time for the joint precedence graph (1) and the joint precedence arc set (2) are found using the following definitions:

t i = q∈Q p q t iq ∀i ∈ I (1) E = q∈Q E q \{reducdant arcs} (2)
Eq. ( 1) creates a weighted average task time for each disassembly task in the final joint precedence graph that considers all EOL states Q. E is the arc set for the final joint precedence graph, and is the union of all arcs in the individual EOL state precedence graphs minus any redundant arcs. Redundant arcs are created when there is a task failure, or zero task time, and a task is skipped because it is no longer required for disassembly.

Method for generating a stochastic joint precedence graph

Disassembly is typically a manual removal process, which makes many of the task times a random variable. We assume that each task time that is a random variable will be normally distributed with a known mean and variance. The normally distributed task times will be referred to as iq with variance 2 iq , for all EOL states q. The weighted task time and variance can be calculated using Eqs.

(3) and (4), respectively.

i = q∈Q p q iq ∀i ∈ I (3) 2 i = q∈Q p q (( iq -i ) 2 + 2 iq ) ∀i ∈ I (4)
Eqs. ( 3) and ( 4) are derived based on a mixture distribution of normal variables. When the task time is zero for a certain EOL state, there will be no variance associated with that particular task time in that state and the variance will very simply be zero. It is obvious that the mixture of two or more normal random variables is not always normally distributed; however, we use the calculation of a variance in (4) to show the relative spread of having a mixture of two or more task times.

Example: joint precedence graph generation

Fig. 1 contains a liaison graph and a cross-section view of a pen taken from a paper by De Fazio and Whitney [START_REF] Fazio | Simplified generation of all mechanical assembly sequences[END_REF]. Based on the information taken from their paper, a theoretical baseline disassembly precedence graph can be created and theoretical disassembly task times added for each component, shown in Fig. 2. For simplicity, the disassembly task for each component, represented as nodes in the precedence graph, will include all tasks to remove that particular component from the core. It is possible that multiple tasks are required for component removal but all task times pertaining to removal of a single component will be summed together into a single time. For this theoretical pen disassembly example, complete disassembly is required. The values shown next to each node in Fig. 2 are the baseline disassembly task times. The body of the pen will be the only remaining component after all other components are removed so the 1 time unit (tu) is the amount of time it will take the worker to place the body in a bin or receptacle. EOL data is collected for pens that will be disassembled and the following is the quality information:

• 25% chance that the cap is missing for disassembly, resulting in a zero task time. • 20% chance that the head is damaged so disassembly time increases from 2 tu to 4 tu. • If the head is damaged then 50% chance the tube is damaged and the disassembly time increases from 1 to 3 tu.

Since there are 3 (x = 3) components that have a disassembly task different from the baseline and there are 2 x = 8 possible precedence graphs to consider. However, since one of the EOL possibilities has a dependency between 2 components (head and tube), the 8 possible precedence graphs decrease to 6. The coded EOL condition table is shown in Table 1(a). If there is a 1 in the cell, then that component has the baseline precedence graph disassembly task time and if there is a 0, then that component has the altered disassembly time. Table 1(b) shows the actual probability of the EOL condition for each component. The cap and head components have probabilities that are all less than 1 because they are independent quality types but the tube has a probability of 0.5 or 1 because if the head is the baseline EOL state, the probability that the tube is in a state other than the baseline is 0. This adjustment allows the total probability of the 6 EOL states to add up to 1 since 2 of the EOL states are eliminated due to the dependency between two components. Table 2 shows the probability of each EOL state and the corresponding disassembly task time for the EOL state. EOL state 1 is the baseline precedence graph and there is a 60% chance the pen will be in this EOL state. This EOL data can be incorporated into individual precedence graphs for each EOL state, shown in Fig. 3. Each precedence graph in the table has a number in parenthesis below that corresponds to its EOL state, and a probability of that EOL state. The cap node for EOL states 4, 5, and 6 are dashed because this is essentially a skipped task that results in 0 task time. After using Eqs. ( 1) and (2) to calculate the weighted average disassembly task time for each component i, t i , and E, a joint precedence graph can be created, shown in Fig. 4. The disassembly task time for the cap decreased and the disassembly task times for the head and tube increased to reflect all the possible EOL states the pen can be in for disassembly. Now let us assume that each task time that is non-zero is a normally distributed variable that has a standard deviation equal to 0.5 tu, or a 0.25 tu variance. Table 3 reflects the change that each non-zero task time has a variance of 0.25 tu. Applying Eqs. ( 3) and (4) from Section 2.2, we can generate a stochastic joint precedence graph, shown in Fig. 5.

The stochastic joint precedence graph in Fig. 5 is very similar the joint precedence graph in Fig. 4, except now each task time is a normal variable with a mean and variance.

Disassembly line balancing algorithm

Now that we can create a disassembly joint precedence graph, we introduce a line balancing algorithm that can use the joint precedence graph as an input. This is an example of one line balancing formulation but many can be used where a single precedence graph is utilized.

Deterministic disassembly line balancing formulation

Notation:

• Tasks: i ε I • Workstations: k ε K • Task Time: t i • Predecessor Set: P(i)
• Tool Sharing Set: P(tool i)

• Maximum Tasks Allowable: Z

Objective function:

min max i∈I X ik t i - i ∈I i k t i ∀k, k ∈ K and k < k (5)
Subject to:

k∈K X ik = 1 ∀i ∈ I (6) X ik ≤ k h=1 X i h ∀i, i ∈ I, i / = i , k ∈ K and i ∈ P(i) (7)
X ik = 0, 1 ∀i ∈ I and k ∈ K (8)
The objective function (5) minimizes the maximum difference in total workstation time for a given number of workstations K. X ik is Fig. 5. Stochastic joint precedence graph for pen disassembly.

the binary decision variable that is 1 if task i is assigned to workstation k and 0 otherwise. Constraint ( 6) is an occurrence constraint where every disassembly task must be assigned a workstation and only one workstation. Constraint [START_REF] Altekin | Disassembly line balancing with limited supply and subassembly availability[END_REF] are the precedence constraints [START_REF] Baybars | A survey of exact algorithms for the simple assembly line balancing problem[END_REF] where the workstation assignment of task i l will be at the same workstation or an earlier workstation of task i, based on the predecessor set P(i). Constraint ( 8) is a non-divisibility constraint that requires task i to not be split amongst multiple stations.

X ik = X i k ∀i, i ∈ P(tool i), k ∈ K (9) i∈I X ik ≤ Z ∀k ∈ K (10)
Constraint ( 9) is a tool sharing constraint where task i and i l are tasks in a set P(tool i) and have a common disassembly tool required to complete the task. In order to reduce cost for additional tooling, the disassembly tasks in this set must be completed at the same workstation. Constraint ( 10) is a total workstation task constraint that limits the number of disassembly tasks Z that can be performed at a single workstation. Since disassembly is typically manual, if workers are overloaded with too many tasks then the quality of their work can decrease and the time to complete tasks can increase.

To model the min-max function in CPLEX, the objective function is transformed and 2 additional constraints are added.

Objective function:

min Y (5a) Subject to: i∈I X ik t i - i ∈I X i k t i ≥ -Y ∀ k, k ∈ K, and k < K (5b) i∈I X ik t i - i ∈I X i k t i ≤ Y ∀ k, k ∈ K, and k < K (5c) (6)
, ( 7), ( 8), ( 9), [START_REF] Altekin | Profit-oriented disassembly-line balancing[END_REF] Objective function (5a) minimizes a float variable Y that is in constraint (5b) and (5c). Constraint (5b) ensures that the difference between total workstation time is greater than or equal to the negative of Y and constraint (5c) ensures that difference between the total workstation time is less than or equal to the positive value of Y.

Stochastic disassembly line balancing formulation

The stochastic disassembly line balancing formulation is very similar to the deterministic line balancing formulation in Section 3.1, except now we will take into consideration that each task time is a normal random variable and the total time variance for each workstation needs to be below some predetermined value. The objective in [START_REF] Tripathi | Real world disassembly modeling and sequencing problem: Optimization by Algorithm of Self-Guided Ants (ASGA)[END_REF] is to minimize the maximum difference between the total workstation times. The only change in nomenclature is now task times, V i , are normal random variables with mean i variance 2 i . Objective function:

min max E i∈I X ik V i - i ∈I X i k V i ∀k, k ∈ K, and k < k (11)
Subject to:(6), ( 7), ( 8), ( 9), [START_REF] Altekin | Profit-oriented disassembly-line balancing[END_REF]

) i X ik 2 i ≤ ∀k ∈ K (12)
Constraint ( 12) is a total workstation variation constraint that requires each workstation to have a total variance below some value ˙. We take advantage that the sum of normally distributed variables is a normal distribution with mean that is the sum of all means and a variance that is the sum of all variances. For constraint [START_REF] Koc | Two exact formulations for disassembly line balancing problems with task precedence diagram construction using an AND/OR graph[END_REF], the value chosen for ˙ will dictate how large the total workstation variation can be. The expectation can be removed from the objective function, based on Proposition 1, and we can transform the objective into two constraints so we can plug our formulation into CPLEX.

Objective function:

min M (11a) Subject to: i∈I X ik i - i ∈I X i k i ≥ -M ∀ k, k ∈ K, and k < K (11b) i∈I X ik i - i ∈I X i k i ≤ M ∀ k, k ∈ K, and k < K (11c) (6)
, ( 7), ( 8), ( 9), ( 10), ( 12)

Proposition 1. The E i∈I X ik V i - i ∈I X i k V i ∀i k, k ∈
K, and k < K, where i and i are the means for normally distributed variables, is equivalent to i∈I X ik ii ∈I X i k i .

Proof. The expected value of a normal distribution is its mean, , and the difference between two normally distributed variables is a normal distribution with the mean being the difference in means. Therefore,

E[ i∈I X ik V i - i ∈I X i k V i ] = E[ i∈I X ik V i ] - E[ i ∈I X i k V i ] = i∈I X ik i - i ∈I X i k i . ᮀ

Disassembly line balancing: laptop example

The following laptop example is more realistic example compared to the pen example because laptops have components and modules that have tangible value for reuse and aftermarket sales. The disassembly times and EOL data are hypothetical and used to show the benefit of line balancing using a joint precedence graph versus a single precedence graph that does not consider every EOL state and component/module reuse possibilities. The laptop example has been used in a couple of papers as a demonstration of assembly and disassembly [START_REF] Hu | Assembly system design and operations for product variety[END_REF][START_REF] Riggs | Disassembly liaison graphs inspired by word clouds[END_REF]. The exploded view of the laptop was found on dell's website. The dell laptop has 13 components, as shown in Fig. 6. A disassembly precedence graph was generated by analyzing the connections between the components and disassembly times estimated. The baseline precedence graph for the dell laptop is shown in Fig. 7.

The numeric values next to each node are the total disassembly task time for each task. An end node is added with 0 time so when generating the different line balances, there is a common end task. The following is the collected EOL quality information:

• 30% chance component G (battery) is missing which results in a 0 task time. • 35% chance component I (hard drive) is missing which results in a 0 task time. • 20% chance component E (system board) has stripped screws so the time to remove E from the computer base increases from 5 tu to 8 tu. • 30% chance components A (display), B (hinge), C (keyboard), and D (palm rest) do not need to be disassembled and can be left as a module to directly be reused in the computer aftermarket.

Reusing the module will result in a 0 task time for A, B, and C but D will still take 4 tu.

The EOL information was derived to somewhat be based on realistic EOL scenarios, although the percentage of each EOL state and the task times are not based on an actual case study. It is very possible that consumers will keep the battery as a backup for their next computer if it is compatible, and the hard drive may be removed to save all data. Additionally, these components may be removed and sold in the aftermarket by individuals. It is assumed the system board is screwed into the computer base and small screws can present difficulties for workers to undue. The modularization of components A, B, C, and D for reuse is to simulate how different types of reuse are possible for the same components. If each individual component is disassembled they can be inspected and determined fit or not for reuse, and if not fit they can be recycled; however, there is a chance the entire module is found fit as is and instead of completely disassembling the components, they can be left together and go straight into reuse.

There are k = 4 different components with various EOL states, so there are 16 possible EOL disassembly precedence graphs. Instead of showing each disassembly precedence graph, the information is shown in Table 4 with all relevant information.

Similar as before with the pen example, EOL state 1 is the baseline case and the probability of this EOL state presenting itself for disassembly is 20.7%, the highest of all EOL states. The data in Table 4 can be used to create a joint disassembly precedence graph for the laptop example, shown in Fig. 8. Versus the baseline case, the disassembly task time in the joint precedence graph for components A, B, C, I, and G decreased, while the disassembly task time for component E increased.

It is very difficult to include all 16 EOL states for the laptop example into a single line balancing algorithm and have it solve in a timely manner. For other realistic disassembly examples, the number of EOL states can increase dramatically as more EOL condition data and reuse possibilities are included into calculating EOL states. To show the effectiveness of creating a disassembly joint precedence graph for line balancing, we use CPLEX to determine two optimal line balances. The first line balance uses the joint precedence graph as an input and the second line balance uses the EOL state that has the highest probability of presenting, which in most cases will be the baseline precedence graph. We then created a simulation for each line balance and compared the average throughput between the two. When running CPLEX to obtain the optimal line balances, the solution time never exceeded 2 s.

For our simulation we assume the first workstation in the serial line cannot be starved and the final workstation cannot be blocked. We used the line balancing setup from Section 3 that minimized the maximum difference in total workstation time. We simulated the 2 line balances using Simul8 software and ran each simulation 5 times with different randomly generated arrival seeds. Each simulation was run for 2400 h to reach a steady state and the disassembly task times were turned into minutes, for example 2 tu equals 2 min. The time unit does not matter as we find the percent difference in throughput for the different simulations, which is unit less. Table 5 contains the line balancing task assignment for each workstation for the baseline and joint cases, as well as the average throughput and percent difference. We balanced for 4 workstations and this remained constant. Table 6 contains the average workstation utilization percentages (waiting, working, and blocking percents) for the simulation runs of each line balance.

Line balancing sensitivity

The results in Table 5 are for a combination of different probabilities of EOL states but it does not give any sort of understanding of the sensitivity between the probabilities of each EOL state versus balancing against the EOL state that has the highest probability of presenting for disassembly. We created and analyzed a series of line balances and simulations for having each of the 4 previously stated EOL component possibilities (battery and hard drive missing, system board damaged, monitor/keyboard module for reuse) where each condition has a 50%, 40%, 30%, 20%, 10%, and 5% probability of occurring. The same setup was performed for these scenarios as was done before with the simulation results shown in Table 7.

Stochastic line balancing sensitivity

It was assumed before that all task times are deterministic but it is possible that each task can have a small amount of variation due to the manual nature of the work. We performed additional line balances and simulations for the original laptop example except now each non-zero task time is a normal random variable with a mean that is equal to the previously used task times and has a standard deviation of 0.5 tu and 1.0 tu. We used a value of ˙ = 14 tu 2 for the total variance in constraint [START_REF] Koc | Two exact formulations for disassembly line balancing problems with task precedence diagram construction using an AND/OR graph[END_REF] for the stochastic line balancing formulation in Section 3.2. Tables 8 and9 show the mean, variance, and standard deviation for each task using Eqs. ( 3) and ( 4) with a task standard deviation of 0.5 and 1.0 tu, respectively, and the simulation results for the line balances are shown in Table 10. The line balancing task assignment in Table 7 is the same as in Table 4 from the deterministic case because the total workstation variance constraint did not require tasks to be reassigned to different workstations. We also experimented with only 1 task being a normal variable, E (system board) in this case, with a standard deviation of 1.0 to see if the same results would occur as the case with all non-zero task times being normal random variables. The results from this set of experimentations are shown in Table 11.

Discussion

The results in Table 5 show that when creating and simulating a serial disassembly line with 4 workstations, a line balance utilizing the joint precedence graph has an average throughput of 9.72% versus a line balance using the baseline precedence graph. The percent different in each table always compares the increase in throughput the joint line balance has when compared to the baseline line balance. Table 6 shows that the joint precedence graph line balance has an overall higher average workstation working percent, 84.23% versus 307.40%. It is interesting to note that for the baseline line balancing case, a single workstation (in this case workstation 4) will be the bottleneck for the disassembly line with a working percent of 100%. The workload for the joint line balance is more balanced with a working percent range of 83.11-76.85%. This conclusion was fairly constant for the line balancing results in the sensitivity study, where the joint line balance had a relatively small working percent range and a single workstation in the baseline line balance is the bottleneck.

The actual task assignment between the baseline and joint line balances have some identical task assignments with a few different task assignments that creates the more balanced work load for the joint line balance. For both the joint and baseline line balances, workstation 1 is assigned tasks G, A, and D; workstation 2 is assigned C and B; workstation 3 is assigned I; and workstation 4 is assigned M, E, and H. The difference in task assignments are workstation 1 for the joint case is also assigned task L while the baseline balance is not; workstation 2 for the joint case is assigned task J and workstation 2 for the joint case is assigned L and K, but not J; workstation 3 for the joint case is assigned K and F while the baseline balance for workstation 3 is assigned task J, but not tasks K and F; and workstation 4 for the baseline case is assigned F and the joint case is not. Due to the work balancing for the joint line balance, each workstation shifts between being the pacing workstation of the line while the baseline line balance is completely paced by the final workstation.

Line balancing sensitivity discussion

The line balancing sensitivity results in Table 7 are fairly expected. When moving from all EOL states having equal probability of occurring (50%) to diminishing percent of occurrence, the joint line balance decreases in throughput advantage versus the baseline case. The only result that does not match with other results is if every EOL state has a 20% chance of occurring. When each EOL state is 30%, the difference in throughput is 7.31%, for 20% case it drops to 2.29%, and 10% case it rises to 3.03%. When investigating why this result happens, the task assignment for workstation 2 for the joint line balance (tasks A, C, and I) has a 10.5% chance of resulting in a 0 total workstation time due to the task assignment. This result is shown in Table 12 and the EOL states that result in a total workstation time of zero is highlighted in gray.

This observation is typical for disassembly where some tasks should be ignored depending on the EOL state. Therefore, the following important constraint should be added to the model (5a)-(11a) to prevent such situations: j∈O X jk < i∈I X ik ∀k ∈ K, ∀O ∈ O [START_REF] Altekin | Task-failure-driven rebalancing of disassembly lines[END_REF] where O is a family of such sets O that O contains all tasks having 0time for one coded EOL state. After the inclusion of this constraint, task J from workstation 3 was switched with task I to prevent the zero work time from happening. When simulating this new line balance, the average percent difference in throughput was 4.38% in favor of the joint line balance. times are continuous and the variation of workstation time can vary with a higher probability as the standard deviation increases, the interaction of workstations along the disassembly line that can cause waiting and blocking have an infinite amount of possibilities. This should, as expected, lead to a diminished average throughput. For the baseline case, the throughput decreases from 4.598 to 4.572 to 4.479 as the standard deviation increases, and the joint case throughput decreases from 5.045 to 4.881 to 4.702. What was not originally expected is how the decreases in average throughput are more drastic for the joint case versus the baseline case. The average throughput difference between each comparison decreases from 9.72% to 6.76% to 4.98%. The reason for this disparity is the baseline line balances have a single pacing workstation versus the joint line balances have pacing workstations that change throughout the simulation. Tables 13 and14 show how the workstation utilization changes for task time with different standard deviations. Because the baseline line balances have a single pacing workstation, this workstation stays working for close to 100% of the time for all three standard deviation cases (0.0, 0.5, and 1.0 tu) while the workstation working utilization for the joint cases decreases steadily. The range of working percent stays relatively small for each of the joint cases, but the average workstation working percentage decreases as the standard deviation values increase.

The situation where every task is a random variable with a common variance is unlikely, but it is very possible that at least one task time is a random variable. Table 11 in Section 4.2 shows the results for having task E have a standard deviation of 1.0 tu and the throughput results do not change very much when compared to having all task times be deterministic.

Conclusions

This paper considers the variation of components end state and adapts an efficient method from mixed model assembly to plan the disassembly process of such products. Disassembly of products at their EOL and efficient disassembly systems are important topics and a method that is able to consider the full spectrum of EOL conditions and treatment possibilities can be a powerful tool for disassembly decision makers. The goal of this paper is to develop and validate a disassembly joint precedence graph creation method that is able to simultaneously consider all possible EOL conditions and states a product can be in. We also wanted to extend our approach to consider the possibility of stochastic task times and develop a method for generating a stochastic disassembly joint precedence graph.

The achievements of our proposed method is: (1) a method that can handle, through preprocessing, many different EOL states and treatment options for components in the same product into a single precedence graph, and (2) a method that can work with any line balancing algorithm where a single precedence graph is used. The impact our method has on disassembly line efficiency and throughput when using a joint precedence graph versus a precedence graph for only one EOL state is highly dependent on the percentage of EOL states, different treatment options, and task times, both deterministic and stochastic. However, any gain in efficiency for disassembly lines is important, no matter how small. Future work can focus on methods to mitigate the decrease in throughput percent difference when task times are continuous. Since disassembly is mostly manual operations, it makes sense that task times will have some sort of variance associated with them. Additional future work should focus on the possibility that the EOL condition and state probabilities can change over time and how can such changes be handled with respect to reconfiguration of the disassembly line.

Fig. 1 .

 1 Fig. 1. (a) Pen liaison graph and (b) pen cross-section [24].

Fig. 2 .

 2 Fig. 2. Pen disassembly precedence graph.

Fig. 3 .

 3 Fig. 3. EOL disassembly precedence graphs for pen.

Fig. 4 .

 4 Fig. 4. Joint precedence graph for pen disassembly.

Fig. 6 .

 6 Fig. 6. Exploded view of dell laptop and parts list [21].

Fig. 7 .

 7 Fig. 7. Baseline disassembly precedence graph of laptop.

Fig. 8 .

 8 Fig. 8. Disassembly joint precedence graph for laptop example.

Table 1

 1 (a) Coded EOL state and (b) disassembly probability of each EOL state.

	EOL state	Cap	Head	Tube
	(a)
	1	1	1	1
	2	1	0	1
	3	1	0	0
	4	0	1	1
	5	0	0	1
	6	0	0	0
	(b)
	1	0.75	0.80	1.00
	2	0.75	0.20	0.50
	3	0.75	0.20	0.50
	4	0.25	0.80	1.00
	5	0.25	0.20	0.50
	6	0.25	0.20	0.50

Table 2

 2 Probability and time data for each EOL state.

	EOL state	Probability	Cap	Head	Tube
	1 (baseline)	0.600	2	2	1
	2	0.075	2	4	1
	3	0.075	2	4	3
	4	0.200	0	2	1
	5	0.025	0	4	1
	6	0.025	0	4	3

Table 3

 3 Probability and time data for each EOL state.

Table 4

 4 Laptop EOL precedence graph information.

	EOL State	Probability	Disassembly task times (tu)
			A	B	C	D	E	F	G	H	I	J	K	L	M
	1 (baseline)	0.207	4	4	3	5	5	3	2	2	6	5	3	1	2
	2	0.089	0	0	0	5	5	3	2	2	6	5	3	1	2
	3	0.111	4	4	3	5	8	3	2	2	6	5	3	1	2
	4	0.048	0	0	0	5	8	3	2	2	6	5	3	1	2
	5	0.089	4	4	3	5	5	3	0	2	6	5	3	1	2
	6	0.038	0	0	0	5	5	3	0	2	6	5	3	1	2
	7	0.048	4	4	3	5	8	3	0	2	6	5	3	1	2
	8	0.020	0	0	0	5	8	3	0	2	6	5	3	1	2
	9	0.111	4	4	3	5	5	3	2	2	0	5	3	1	2
	10	0.048	0	0	0	5	5	3	2	2	0	5	3	1	2
	11	0.060	4	4	3	5	8	3	2	2	0	5	3	1	2
	12	0.026	0	0	0	5	8	3	2	2	0	5	3	1	2
	13	0.048	4	4	3	5	5	3	0	2	0	5	3	1	2
	14	0.020	0	0	0	5	5	3	0	2	0	5	3	1	2
	15	0.026	4	4	3	5	8	3	0	2	0	5	3	1	2
	16	0.011	0	0	0	5	8	3	0	2	0	5	3	1	2

Table 5

 5 Results of laptop line balancing comparison.

	EOL State	Probability	Line balance				Average throughput	Percent difference
			Wkst 1	Wkst 2	Wkst 3	Wkst 4
	Baseline	20.70%	GAD	CBLK	JI	FMEH	4.598	9.72
	Joint	-	GADL	CBJ	KFI	MEH	5.045

Table 6

 6 Workstation utilization results from simulation.

	Baseline case

Table 7

 7 Results of laptop line balancing sensitivity.

	Component probability	EOL state	Probability	Line balance				Average throughput	Percent difference
				Wkst 1	Wkst 2	Wkst 3	Wkst 4
	50%	Baseline	6.25%	GAD	CBLK	JI	FMEH	4.444	13.07
		Joint	-	GDLM	ABJ	KFI	CEH	5.025
	40%	Baseline	12.96%	GAD	CBLK	JI	FMEH	4.546	12.30
		Joint	-	GADL	BJM	KFI	CEH	5.105
	30%	Baseline	24.01%	GAD	CBLK	JI	FMEH	4.65	7.31
		Joint	-	GADL	BFI	KJM	CEH	4.99
	20%	Baseline	40.96%	GAD	CBLK	JI	FMEH	4.762	2.29
		Joint	-	GDLF	ACI	BJM	KEH	4.871
	10%	Baseline	65.61%	GAD	CBLK	JI	FMEH	4.877	3.03
		Joint	-	GDLF	AMI	CKJ	BEH	5.025
	5%	Baseline	81.45%	GAD	CBLK	JI	FMEH	4.938	1.58
		Joint	-	GDLK	AMI	CFJ	BEH	5.016

Table 8

 8 Data input for the stochastic line balancing formulation, 0.5 tu standard deviation.

	Task	Mean time	Variance	Standard deviation
	A	2.80	3.54	1.88
	B	2.80	3.54	1.88
	C	2.10	2.07	1.44
	D	5.00	0.25	0.50
	E	6.05	2.30	1.52
	F	3.00	0.25	0.50
	G	1.40	1.02	1.01
	H	2.00	0.25	0.50
	I	3.90	8.35	2.89
	J	5.00	0.25	0.50
	K	3.00	0.25	0.50
	L	1.00	0.25	0.50
	M	2.00	0.25	0.50

Table 9

 9 Data input for the stochastic line balancing formulation, 1.0 tu standard deviation.

	Task	Mean time	Variance	Standard deviation
	A	2.80	4.06	2.01
	B	2.80	4.06	2.01
	C	2.10	2.59	1.61
	D	5.00	1.00	1.00
	E	6.05	3.05	1.75
	F	3.00	1.00	1.00
	G	1.40	1.54	1.24
	H	2.00	1.00	1.00
	I	3.90	8.84	2.97
	J	5.00	1.00	1.00
	K	3.00	1.00	1.00
	L	1.00	1.00	1.00
	M	2.00	1.00	1.00

Table 10

 10 Results of laptop line balancing for random variable task times.

	Component standard deviation	EOL state	Probability	Line balance				Average throughput	Percent difference
				Wkst 1	Wkst 2	Wkst 3	Wkst 4
	0.5	Baseline	20.70%	GAD	CBLK	JI	FMEH	4.572	6.76
		Joint	-	GADL	CBJ	KFI	MEH	4.881
	1.0	Baseline	20.70%	GAD	CBLK	JI	FMEH	4.479	4.98
		Joint	-	GADL	CBJ	KFI	MEH	4.702

Table 11

 11 Results of laptop line balancing for only E as random variable.

	Component standard deviation	EOL state	Probability	Line Balance				Average throughput	Percent difference
				Wkst 1	Wkst 2	Wkst 3	Wkst 4
	E = 1.0 only	Baseline	20.70%	GAD	CBLK	JI	FMEH	4.587	9.16
		Joint	-	GADL	CBJ	KFI	MEH	5.007

Table 12

 12 Workstation #2 for 20% EOL state.

	Attribute	Probability	A	C	I	Total work time
	1	0.2070	4	3	6	13
	2	0.0887	0	0	6	6
	3	0.1115	4	3	6	13
	4	0.0478	0	0	6	6
	5	0.0887	4	3	6	13
	6	0.0380	0	0	6	6
	7	0.0478	4	3	6	13
	8	0.0205	0	0	6	6
	9	0.1115	4	3	0	7
	10	0.0478	0	0	0	0
	11	0.0600	4	3	0	7
	12	0.0257	0	0	0	0
	13	0.0478	4	3	0	7
	14	0.0205	0	0	0	0
	15	0.0257	4	3	0	7
	16	0.0110	0	0	0	0

Fig. 9. Comparison of baseline and joint line balances average throughput.

Stochastic line balancing sensitivity discussion

When comparing the task assignments in Tables 5 and10, the baseline line balances are the same in each and the joint line balances are also the same. Table 5 contains results where each task time is deterministic, while Table 10 has task times that are normally distributed with a mean equal to the deterministic times and standard deviations of 0.5 and 1.0 tu, respectively. The percent difference in throughput between the joint versus baseline line balances decreases as the standard deviation increase from 0 (the deterministic case) to 0.5 and 1.0 tu. These results are shown in Fig. 9. These results are intuitive because if the total workstation


Madjid Haddad
email: madjid.haddad@isae.fr

Yves Gourinat

Miguel Charlotte

Equivalence Theory Applied to Anisotropic Thin Plates

Keywords: Equivalence Theory, Thin Plates, Anisotropic Plates, Dynamics

. We provide a numerical analysis of the convergence properties of ET approach while comparing its numerical results with those pre-dicted by the analytical theory and the Finite Elements Method for thin plates.

Introduction

Around the 70s some researchers [START_REF] Absi | La Theorie des Equivalences et Son Application a l'Etude des Ouvrages d'Art[END_REF][START_REF] Vegas | Application de la Theorie des Equivalences a l'Etude d'Une Dalle Biaise[END_REF] were interested in the problem of modelling a slab as a beam lattice or truss. Such a substitution relies on some common engineering viewpoint according which if the slab can be subdivided into slices parallel to its boundaries (like AB or CD on Figures 1 and2), then each slice may be assimilated to a beam with related material and geometrical properties. This coarse structure approach which is suitable for pre-sizing was also assumed to be applicable to the calculation of arch dams (Figure 3) sliced into arcs (AB) and consoles (CD).

Such a practical problem gave birth to the Equivalence Theory (ET) developed by Absi [START_REF] Absi | La Theorie des Equivalences et Son Application a l'Etude des Ouvrages d'Art[END_REF] who specified the general conditions of equivalence in statics for isotropic materials. In statics, the standard equivalence criterion is expressed between the strain energies. That theory assumes that any real body modelling may be substituted by another one that, even though it may possibly have material constitutive laws and geometric properties with no physical meanings (like negative cross sections or Young modulus), is intended to be more advantageous for calculus. In the literature [START_REF] Absi | La Theorie des Equivalences et Son Application a l'Etude des Ouvrages d'Art[END_REF][START_REF] Haddad | Application de la Methode des Equivalences en Dynamique[END_REF][START_REF] Cucchi | Elastic-Static Analysis of Shear Wall/Slab-Frame Systems Using the Framework Method[END_REF], only few investigations have tried to check this assumption. These equivalence analyses were formulated for the statics of different cases of isotropic slabs and the comparisons made with the analytical theory and the Finite Element Method (FEM) analysis have led to very encouraging conclusions. However, the ET has been abandoned in favour of the FEM because the latter is more flexible to deal with structures with arbitrary geometry. This work aims to extend the ET approach to the analysis of anisotropic plates subjected to transverse loads, in statics, and to study the frequencies and the modes shapes in dynamics. To check the validity of the extended formulation, numerical comparisons with analytical solutions in statics and to FEM results in dynamics are made. As mentioned previously, the possible lack of meaning of the resulting equivalent model makes the treatment of such a method on commercial FEM software (Patran/Nastran, Catia...) difficult, if not impossible. Therefore in order to compare our new improvement of the theory and overcome these difficulties, we were led to develop a finite element computer codes in Matlab.

Formulation of the Method in Statics

The Equivalence Theory (ET) aims at replacing a given real body by another arbitrary chosen, and possibly fictitious, one [START_REF] Absi | La Theorie des Equivalences et Son Application a l'Etude des Ouvrages d'Art[END_REF]. The standard equivalence criterion is expressed between the strain energies. Here we will choose the substituting body as a lattice structure with beams elements. We will first express the energy of tensile and bending of an anisotropic plate. Then, according with the ET, we will identify the cross sections, the quadratic moments and the central moments of inertia of beams, by making a comparison between the Elementary Representative Cells (ERC) of the two models.

Expressing Equivalence in Traction

Elastic Energy of an Orthotropic Plate: [4]

Consider the case of a thin, linearly elastic and anisotropic plate, with three planes of symmetry and a constant thickness h. The plate stress strain relations are written in this case as follows: ;

; 2

x x x y y y y x xy xy E E E E G                  (1) with x l E kE   ; y t E kE   ; tl l lt t E k E k E      ; lt G G  Hereabove , ,
x y xy    denote respectively the normal stresses with respect to x and y axes, and the shear stress, while 1 ;

; 2 Suppose now that the plate is subjected to a uniform bi-axial tensile along the x and y axis. Hence its elastic energy reads like

x y x y u v u v x x x y                      (
  2 2 2 2 2 2 p t xx y y xy x y A U E E E G              (3)
where A p represents the ERC surface area of the plate while noticing that the strains are independent on x, y or z. Note that we recover the traction energy of a materially isotropic plates in plane stress for t l E E E   (the same Young's modulus) and

tl lt      (same Poisson's ratio).

Beam Traction Energy

The traction energie of a beam "ij", Figure 4 is:

ij ij W e e e             (4)
where

ij ij e ES l  
is a characteristic traction parameter, while S ij , l e and E are respectively the cross section, length and Young's modulus of the beam.

Expressing the Traction Equivalence

The equivalence is expressed in an elementary representative cell ERC of the two structures.

Consider the cell represented in Figure 5.

The strain energy expressions of the different beams constituting this cell are: Write the equivalence between the two cells, we have then:

   
                       (5)
          2
                                      (6)
Taking AD BD    , and comparing the two terms of the Equation ( 6) we can write:

      2 2 2 4 11 11 11 4 cos 2 cos 2 p x AB AD p x AB AD A E e e e A E             (7)       2 2 2 4 22 22 22 4 sin 2 sin 2 p y AC AD p y AC AD A E e e e A E             (8)     2 2 2 2 12 12 2 2 4 sin cos 1 4 sin cos p A D p AD A G e e A G         (9)
  Replacing AD  in the Equations ( 7) and (8), we obtain:

  2 1 cot 2 p l AB tl A E k      (11)   2 1 tan 2 p t AC lt A E k      (12) 2 2
2 sin cos

p l tl AD A E      ( 13
)
We know that:

  AB AB EA a   ,   AC AC EA b   and   AD AD EA l
  ; where A ij , E ij are respectively the cross section and the Young's modulus of the beam "ij" and a, b, l: lengths of the considered beam, as shown in Figure 5.

Let's consider that the beams which are parallel to the x axis and the diagonal ones have the same Young's modulus as the longitudinal one of the plate, and those which are parallel to the y axis, have the same Young's modulus as the transverse one of the plate.

Finally, we obtain:

  2 1 cot 2 AB lt bh A k     ;   2 1 tan 2 AC lt ah A k     ; 2 2 2 sin cos lt AD abh A k l           .

Bending Equivalence

Plate Bending Energy

In our study we consider the Kirchoff-Love's plates, which suppose that the straight linear elements that are perpendicular to the plate's mean surface still remain so even after deformation. We have:

2 2 2 2 2
; ; .

x y x y w w w z z z x y x y                 (14)
Then the bending and tensional moments are given as follows:

/2 2 2 1 2 2 / 2 / 2
                                        (15)
where

3 12 x x E h D   ; 3 12 y y E h D   ; 3 1 12 E h D   ; 3 12 xy Gh D  .
The bending energy of the plate is given [START_REF] Timochenko | Théorie des Plaques et Coques[END_REF]:

2 2 1 d d d 2 x y x y w w w V M M M x y x y x y                                  , ( 16
)
which gives:

2 2 2 2 2 1 2 2 1 d 2 2 dd xy w w V D D x y x y w w w D D xy x y x y                                                (17)
with:

3 12 l x kE h D  ; 3 12 t y kE h D  ; 3 1 12 tl l k E h D   ; 3 12 lt xy G h D  .

Bending and Torsional Energie of a Beam

The bending and tensional energy contributions of a beam 'ij'(see [START_REF] Absi | La Theorie des Equivalences et Son Application a l'Etude des Ouvrages d'Art[END_REF][START_REF] Timochenko | Théorie des Plaques et Coques[END_REF]) are respectively given as

2 2 2 2 2 2 2 2 1 cos sin 2 cos sin f i j w w W x y w x y                       (18)
where

  ij e ij EIl  
is defined with the quadratic moment I of the cross section S ij and

2 2 2 2 2 2 1 1 sin 2 cos 2 2 2 t w w w W x y y x                            (19)
where e GJl   .

Expressing the Bending Equivalence

Here we will assume that all the beams work in bending and only the longitudinal and diagonal ones work in torsion (Figure 6). This supposition has no effect in the theory itself, since the unique condition to satisfy is to conserve the total strain energy. We can write then:

2 2 2 2 2 1 1 2 2 fAB fCD AB AB w w W W x y x                        2 2 2 2 2 1 1 2 2 fAC fBS AC AC w w W W x y y                        2 2 2 2 2 2 2 2 1 cos sin cos sin 2 fAD AD w w w W x y x y                       2 2 2 2 2 2 2 2 2 1 cos sin cos sin 2 BC w w w W x y xy                     
The equivalence is given by equating the two strain energies:

2 2 2 2 fAB fBC fAC fAD tAB tAC W SV W W W W W W        (20)
                                                                                                     2 2 2 2 2
sin cos

AB AC w x y w w x y x y                                 (21) Let now AD BD     
we can write:

2 2 2 2 2 2 4 2 2 2 4 cos 2 cos 2 x AB AD x AB AD SD w w w x x x SD                                       (
                                      (23) 2 2 2 2 2 2 1 2 2 2 2 1 2 2 2 sin cos 1 2 sin cos AD AD w w w w SD x y x y SD                   (24)
                                                      (25)
Replace the value of AD  in the different equations, we have then: Here I ij and J ij denote respectively the second and polar moments of the cross section area of the beam "ij".

2 1 cot 2 2 x AB SD SD      ; 2 1 tan 2 2 y AC SD SD      ; 1 2 2 1 2 sin cos AD CB SD         ;
By substituting these values in (36), it comes then: The theoretical transverse displacement solution of this orthotropic plate problem at any given point with coordinates (x, y) is given in [START_REF] Abrate | Inpact in Composite Structures[END_REF] page 59 as m n mn    (28).

  3 2 1 cot 24 AB tl bh I k     ;   3 2 1 tan 24 AC lt ah I k     ; 3 2 2 24 sin cos tl AD abh I k l           ;   3 1 6 tl l AB AC lt E abh J J k a b G            .

Static Validation

  1 1 π π , s i n s i n mn n m m x n y w x y W L l                    (
mn mn q W D m D D m n r D n r a               (27)
The deformed shape obtained after loading is given in Figure 8.

A comparison between the theoretical deformation and the one given by the ET, along the two median lines, X = 0.5 and Y = 0.4, is shown in Figure 9.

To evaluate the accuracy of this method, a comparison with the FEM software Patran/Nastran is done for the same size of the mesh and number of nodes as in ET. Quadratic elements are chosen to mesh the thin plate structure in Patran/Nastran. The different results are resumed in Table 1, and compared to the theoretical displacement in the plate center w theoric = 0.0142368 m. While a better accuracy is reached with the FEM method, as expected, very good agreements are observed between the FEM and ET and even the accuracy with the later still remains good for relatively coarse mesh sizes. Moreover, the convergence to the theoretical solution is also observed with increasing the mesh density.

Simply Supported Plate Submitted to Concentrated Load Applied in Its Center

Let us consider the previous thin plate is now loaded at its center by a concentrated force of magnitude P = 0.01 N.

The theoretical expression of the displacement solution given in ( 26) and ( 27) stay the same, except that 0 0 4 sin π sin π

mn x y P q m n Ll L l              ( 29
)
The deformation caused by this load is represented in Figure 10. The related deformation obtained by the ET approach is represented on Figure 11 and a comparison with the theoretical solution is performed along the two median lines in Figure 11.

The Table 2 provides the related numerical values with in addition the results obtained by the FEM theory obtained for the same size mesh, (same number of nodes). The calculation of the theoretical displacement in the plate center gave w theoric = 0.00132144758109 m. Once again, we observe that the two methods converge to the theoretical solution, when we increase the mesh size. Unlike the first example, we observe here (as the main result of the statics section) that the ET provides better results than those of the FEM.

Application of the Equivalence Theory in Dynamics

In statics we have observed that the ET provides very good results, which are in agreement with the literature [1-3], and we have obtained the same conclusions with anisotropic plates. Now we will test this method in dynamics.

In our study we consider the mass conservation of the system. By equating the overall masses of the continuous plate and lattice models, and supposing that a homogeneous mass density in the lattice structure, we have: To calculate the inertial moment of the beams we have to precise a shape for the cross sections. We choose circular ones, having the section from the static equivalence; we calculate first the beam cross section radii π The plate considered here is the same as in the first example, we have just clamped it in its four sides.

In dynamics we preferred to compare the ET results to the FEM ones. We have choose a mesh of 200 × 160 elements in Nastran/Patran software, and we defined the plate as a shell element, these was to approximate better the real solutions.

In what follows we will show the first hundred "100" frequency values obtained for different mesh sizes, and for two representations of the mass matrix, the first is the consistent mass matrix, which is the same expressed in the classic "FEM" formulation, and the second one is the lumped mass matrix, whish considered that the mass and the inertial moment is equitably distributed in the nodes of the beams. Moreover we will have for the translation the half mass working in the three direction "u, v, w" and the half inertial moment working in the three rotating directions " , ,

x y z    " .Note that the bending isn't considered here, this approach is the same as the Patran/Nastran one for the lumped mass matrix.

Figure 12 shows the first hundred eigenvalues of the plate obtained by the ET and confronts them with the values obtained by finite element method. We see very clearly that when the mesh size in the ET is increasing, the curve obtained by this method tends to approximate that obtained by FEM. We also note that even with a coarse size mesh, the first eigenvalues are very close to those obtained by the FEM. The evolution of the relative discrepancy between the ET results versus the FEM ones is plotted in Figure 13. In fact, we observe that the first five "5" frequencies obtained are very close to those obtained by FEM even when we use small size mesh. The convergence as the number of elements is increasing is very remarkable, we note that the first hundred "100" frequencies are obtained with a maximum error of 25% for a mesh of (20 × 16) to 12% for a mesh of (30 × 24). The convergence is also observed in the case a distributed (lumped) mass matrix. Indeed we notice that with the mesh size of (20 × 16) we have a maximum error of 86%, this error is reduced to 16% for a mesh size of (30 × 24). We clearly note that the best results are obtained for a consistent weight distribution.

Mode Shapes

Figure 14 illustrates the first ten mode shapes obtained by the ET approach. These mode shapes are exactly the same than those found in literature [START_REF] Leissa | Vibrations of Plates[END_REF].

Conclusions

Heretofore the equivalence theory has been expressed only for isotropic static problems. The formulation is based in equating the strain energies of two systems. The researchers [START_REF] Absi | La Theorie des Equivalences et Son Application a l'Etude des Ouvrages d'Art[END_REF][START_REF] Vegas | Application de la Theorie des Equivalences a l'Etude d'Une Dalle Biaise[END_REF][START_REF] Cucchi | Elastic-Static Analysis of Shear Wall/Slab-Frame Systems Using the Framework Method[END_REF] were interested by the calculation of deformations, bending moments, convergence with increasing mesh size. They found that in statics this method is very robust, it gives very good results with small mesh sizes and the results are converging when we increase the mesh size.

This method may provide a sheap, sufficient reliable, and convient approach to treat complex structural systems as the computer storage requirement and running times are small compared with those of other numerical methods.

In order to demonstrate the validity of the procedure for anisotropic static problems, two example problems, relevant to plane stress (thin plates) were examined by using the ET approach and the FEM with different grid arrangement. The results are in good agreement with the analytical solutions even in the case of a small number of nodes (coarse mesh). Most importantly, we have also noted that these results are very near and even in some cases better than those predicted by the FEM.

To extend this procedure to dynamics, we have just formulated in addition a simple mass conservation rule, and choose a circular form for our beams. In order to validate this procedure, we compared the results to those of the FEM with a relatively fine mesh. We have not made a comparison with analytical solutions because this latter proposes approximate solutions concerning only the bending of the plates, while the FEM gives the different vibration modes (bending, torsion, traction…). We considered the problem of a thin clamped plate, and found that the results are in good agreement with the FEM solutions. We note also that the formulation of the mass as a consistent one gave best results; this results are explained by the fact that the mass is well distributed in the ERC and approaches the homogeneous distribution in the real cell. The lumped representation of the mass is found to be convenient in the high computational speed, indeed the mass matrixes are diagonals, and so easier to invert, moreover we observe that its results remain close to those of a consistent formulation. Thus the user can be free to choose between very accurate results and high computational speed.

A future work could be to deal with 3D problems. We have created the connectivity table of different 3D "ERC", and begin the solving of the problem concerning the expression of the stiffness matrix in 3D. The analysis performed for the 2D models can be extended to the 3D ones and also addresses other problems like fatigue and blucking.

Figure 1 .

 1 Figure 1. Rectangular slab subdivided into slices parallel to its boundaries.

Figure 2 .

 2 Figure 2. Lozenge slab subdivided into slices parallel to its boundaries.

Figure 3 .

 3 Figure 3. Arch dam subdivided into arcs.

2 )

 2 are respectively the related normal and shear strains; besides E l and E t : are respectively the longitudinal and transversal Young's moduli, lt  and and G lt is the Coulomb's coefficient.

Figure 4 .

 4 Figure 4. A beam representation in the plane.

Figure 5 .

 5 Figure 5. ERC "square with diagonals".

Figure 6 .

 6 Figure 6. ERC «square with diagonals».   1

3. 1 . 7 ),

 17 Simply Supported Plate Submitted to Uniform Distributed Load (Figure Consider that the plate has the following geometrical and material properties: length L = 1 m; width l = 0.8 m, and thickness h = 1 mm, and transversal Young's moduli: E l = 11.10 9 Pa and E t = 11.10 6 Coulomb's coefficient: G lt = 40.10 6 Pa. That plate is submitted to a uniform constant load q 0 = 1 Pa.

Figure 7 .

 7 Figure 7. Simply supported plate submitted to uniform distributed load. r L l  and

Figure 8 .Figure 9 .

 89 Figure 8. Deformation under a uniform distributed load.

  respectively the cross section and length of the e th beam, N is the number of the beams,  the mass density of the beams, p  the mass density of the plate, L, l and h are respectively the length, width and thickness of the plate.

Figure 10 .

 10 Figure 10. Deformed shape under a before concentred centred load.

Figure 11 .

 11 Figure 11. Superposition of the results among the median lines (A-X = 0.5 ; B-Y = 0.4).

Figure 12 .Figure 13 .

 1213 Figure 12. Comparison of the first hundred frequencies obtained by the ET and those provided by the software Patran/ Nastran. (a) Lumped mass matrix; (b) Consistent mass matrix.

Figure 14 .

 14 Figure 14. Mode shapes of orthotropic clumped plate.

Table 1 . Comparison of the central deformation for a sim- ply supported plate under a uniform loading.

 1

	mesh size	FEM deformation [m]	EM deformation [m]	Error FEM (%)	Error ET (%)
	10 × 8	1.405E-02	1.340E-02	-1.34%	-4.59%
	20 × 16	1.422E-02	1.378E-02	-0.12%	-3.11%
	30 × 24	1.425E-02	1.390E-02	0.08%	-2.46%

Table 2 . Comparison of the central deformation for a sim- ply supported plate under a concentrated load applied in its center.

 2

	mesh size	FEM deformation [m]	EM deformation [m]	Error FEM (%)	Error ET (%)
	10 × 8	1.417E-03	1.291E-03	9.78%	-2.34%
	20 × 16	1.367E-03	1.321E-03	3.49%	-0.02%
	30 × 24	1.348E-03	1.325E-03	1.77%	0.24%