JP6446740B2 - Simulation apparatus, simulation method, and program - Google Patents

Description

  The present invention relates to a simulation apparatus, a simulation method, and a program.

It is known to use a yield function such as Yld2000-2d in the simulation of molding by the finite element method (see, for example, Patent Document 1). The yield function includes material property parameters set by the material properties. For example, Yld2000-2d includes a material characteristic parameter αi (i = 1 to 8) and an order M.
Japanese Patent Application Laid-Open No. 2010-61249

  However, conventionally, since a set of material characteristic parameters is used for one material, it cannot be said that anisotropic hardening is sufficiently taken into account, and a deviation from the experiment has occurred.

  According to a first aspect of the present invention, there is provided a simulation apparatus that executes a forming simulation for an object of a predetermined material using a predetermined yield function including a material characteristic parameter, wherein the material has a reference plasticity. A parameter acquisition unit that acquires the relationship between the strain value and the material property parameter value that corresponds to the reference plastic strain, and the material property parameter value that corresponds to the reference plastic strain value when a virtual displacement is applied to the object Are sequentially used as a yield function, and a deformation calculation unit that calculates the deformation state of the object.

  According to a second aspect of the present invention, there is provided a simulation method for executing a forming simulation on an object of a predetermined material using a predetermined yield function including a material characteristic parameter, wherein The parameter acquisition stage for acquiring the relationship between the strain value and the material property parameter value corresponding to the reference plastic strain, and the material property parameter value corresponding to the reference plastic strain value when a virtual displacement is applied to the object Are sequentially used as a yield function to calculate a deformation state of the object.

  According to a second aspect of the present invention, there is provided a program for controlling a computer that executes a molding simulation for an object of a predetermined material using a predetermined yield function including a material characteristic parameter, The parameter acquisition stage for acquiring the relationship between the value of the standard plastic strain and the value of the material property parameter corresponding to the standard plastic strain, and the material property corresponding to the value of the standard plastic strain when a virtual displacement is given to the object The computer is caused to execute a deformation calculation step of calculating the deformation state of the object by sequentially using the parameter values as the yield function.

  The summary of the invention does not enumerate all the features of the present invention. In addition, a sub-combination of these feature groups can also be an invention.

1 schematically shows a test system used in the present embodiment. FIG. 2 is a schematic view showing a state where a circular tube test piece 10 is deformed by a biaxial bulge test apparatus 100. The functional block of the simulation apparatus 200 is shown. An example of the flowchart of operation | movement of the simulation apparatus 200 is shown. It is a flowchart which shows the detail of step S30 of FIG. The stresses σ _x and σ _y for each standard plastic strain ε ₀ ^p are shown. A set of values of the discrete material property parameter αi corresponding to the value of the discrete reference plastic strain ε ₀ ^p is shown. The parameters of each identified function are given. It is a flowchart which shows the detail of step S40 of FIG. The hydraulic bulge forming used in Example 2 is schematically shown. Each material characteristic parameter αi used as a comparative example is shown. The comparison between the experimental value of the true plate thickness strain at the apex | ε _z | versus the internal pressure P diagram and the calculated value by the finite element method is shown. It is thickness distortion distribution along the meridian of 0 degree, 45 degrees, and 90 degrees directions from a rolling direction. It is thickness distortion distribution along the meridian of 0 degree, 45 degrees, and 90 degrees directions from a rolling direction. It is thickness distortion distribution along the meridian of 0 degree, 45 degrees, and 90 degrees directions from a rolling direction.

  Hereinafter, the present invention will be described through embodiments of the invention, but the following embodiments do not limit the invention according to the claims. In addition, not all the combinations of features described in the embodiments are essential for the solving means of the invention.

  FIG. 1 schematically shows a test system used in this embodiment. The test system includes a biaxial bulge test apparatus 100 and a simulation apparatus 200. The test system identifies a material model in consideration of anisotropic hardening for the material using the result of the biaxial bulge test apparatus 100 testing the material of interest, and the simulation apparatus 200 uses the material model. The purpose of this is to obtain a molding simulation result in consideration of anisotropic curability.

  The biaxial bulge test apparatus 100 independently applies an axial force T in the axial direction φ and an internal pressure P in the radial direction θ to the circular tube test piece 10. The biaxial bulge test apparatus 100 includes a main body 102, a pair of drive units 120 and 122, a pressure adjustment unit 130, a pair of cameras 140 and 142, and a controller 104 that controls them.

  Each of the pair of drive units 120 and 122 includes gripping portions 110 and 112 that grip the axial end portion of the circular tube test piece 10. With the gripping portions 110 and 112 gripping the axial end of the circular tube test piece 10, the pair of driving portions 120 and 122 moves in the axial direction φ, so that the axial force T as tension or compression force is changed to the circular tube. The test piece 10 is given.

  The pressure adjusting unit 130 is connected to the hollow portion 12 of the circular tube test piece 10 through a flow path 132 provided in the driving unit 122. The pressure adjusting unit 130 contains a fluid such as oil. In a state where the pressure adjusting unit 130 is connected to the hollow portion 12 of the circular tube test piece 10 through the flow path 132, the fluid is pressurized to apply an internal pressure P in the radial direction θ to the inner wall of the hollow portion 12.

  The pair of cameras 140 and 142 images a predetermined region on the surface of the circular tube test piece 10 from different directions. In this case, a pattern in which black ink is sprayed on a white background is formed in the region of the tube test piece 10. By imaging this pattern with the cameras 140 and 142 and observing changes in the positions of the respective black spots included in the pattern, how much the region of the tube test piece 10 is deformed in the axial direction φ and the radial direction θ. It can be detected. In this case, since the pair of cameras 140 capture images from different directions, the deformation can be detected three-dimensionally based on the principle of stereoscopic vision.

FIG. 2 is a schematic view showing a state where the circular tube test piece 10 is deformed by the biaxial bulge test apparatus 100. The stress in the axial direction φ generated in the tube test piece 10 is expressed as σ _φ , the strain as ε _φ , the stress in the radial direction θ as σ _θ , and the strain as ε _θ .

Although the axial force T and the internal pressure P can be controlled independently, the circular tube test piece 10 is deformed momentarily in the axial direction φ and the radial direction θ. From T and internal pressure P, stresses σ _φ , σ _θ , strains ε _φ , ε _θ are calculated. The controller 104 calculates the stresses σ _φ , σ _θ , strains ε _φ , ε _θ and feeds back the axial force T and the internal pressure P for servo control.

Controller 104, the ratio of the stress sigma _phi: sigma while maintaining _theta by increasing these, stress sigma _phi on reaching each strain plurality of reference plasticity predetermined epsilon ₀ ^p, the sigma _theta, reference plastic strain The stresses σ _x and σ _y are output to the simulation apparatus 200 in association with ε ₀ ^p . The reference plastic strain ε ₀ ^p will be described later. The controller 104 also outputs to the simulation apparatus 200 the result of a similar test performed by changing the stress ratio σ _φ : σ _θ .

  A yield function is used to represent the elastoplastic deformation behavior of a material. Various yield functions have been proposed, and each yield function includes a material characteristic parameter set in accordance with the material. The inventors have found that when the material property parameters included in the yield function are expressed as a function of the standard plastic strain, the anisotropic hardening behavior can be formulated by the yield function so as to match the experiment. Based on the discovery, material characteristic parameters are identified, and a simulation apparatus 200 that simulates deformation of a test material of the material using the identified material characteristic parameters has been found.

  FIG. 3 shows functional blocks of the simulation apparatus 200, and FIG. 4 shows an example of a flowchart of the operation of the simulation apparatus 200. An example of the simulation apparatus 200 is a personal computer.

  The simulation apparatus 200 includes a stress value acquisition unit 220, a parameter identification unit 230, a deformation calculation unit 240, and a storage unit 250. The stress value acquisition unit 220 and the parameter identification unit 230 constitute a parameter acquisition unit 210.

The stress value acquisition unit 220 obtains the biaxial stress values σ _x and σ _y when the biaxial tensile test is performed so that the biaxial stress ratio is different from the controller 104 of the biaxial bulge test apparatus 100 in advance. A plurality of reference plastic strain values ε ₀ ^{p determined} in association with each other are acquired (S20). The parameter identification unit 230 identifies the value of the material characteristic parameter for each of a plurality of reference plastic strain values based on the plurality of biaxial stress values σ _x and σ _y (S30).

The deformation calculation unit 240 executes a forming simulation by the finite element method by sequentially using the material property parameter αi corresponding to the reference plastic strain value ε ₀ ^p (S40). The solution of the finite element method may be any solution method such as a static implicit solution method or a static explicit solution method.

The controller 104 and the simulation apparatus 200 may be integrated. Further, each function of the simulation apparatus 200 may be executed by installing a software program stored in a medium such as the CD-ROM 202. Further, the stress value acquisition unit 220 may acquire the values of the biaxial stress values σ _x and σ _y and the reference plastic strain value ε ₀ ^p from a database or another PC. Further, the parameter identification unit 230 may be incorporated in a device different from the simulation device 200. In this case, the parameter acquisition unit 210 may acquire the material characteristic parameter from the other device and pass it to the deformation calculation unit 240.

  FIG. 5 is a flowchart showing details of step S30 in FIG. First, a target yield function is selected (S31). For example, the yield function may be set by a keyboard operation from the user, or candidates are stored in advance in the storage unit 250, the candidates are displayed on the display of the simulation apparatus 200, and selected by a keyboard operation from the user. May be. Alternatively, any candidate may be preferentially displayed or automatically selected according to the type of material.

  The yield function is preferably one that can describe anisotropy, and examples thereof include Yld2000-2d, Hill's quadratic yield function, Goto's quartic yield function, and Hosford's yield function. Each of these yield functions includes an unknown material characteristic parameter (hereinafter sometimes referred to as αi) whose numerical value is to be identified by the material.

The parameter identification unit 230 uses the value of the reference plastic strain ε ₀ ^p acquired by the stress value acquisition unit 220 and the values of stress σ _x , σ _y to use the stress points (σ _x , σ _{y) of the} isoplastic work surface. ₎ To identify the material characteristic parameter αi so that the yield surface for each value of the standard plastic strain ε ₀ ^p matches (S32). In this case, the parameter identification unit 230 uses a mathematical method such as a least square method, a material mechanical relational expression such as a dimensionless uniaxial tensile stress such as 0 ° or 90 ° from the rolling direction, or a combination thereof. Use. Thereby, the parameter identification unit 230 identifies a set of values of the discrete material property parameter αi corresponding to the value of the discrete reference plastic strain ε ₀ ^p , respectively.

The parameter identification unit 230 further selects a functional form of the material property parameter αi for the continuous reference plastic strain ε ₀ ^p (S33). For example, the function form is stored in advance in the storage unit 250, the candidate is displayed on the display of the simulation apparatus 200, and is selected by a keyboard operation from the user. Alternatively, any candidate may be preferentially displayed or automatically selected according to the yield function or the like. The function form includes parameters to be identified according to the material.

The parameter identification unit 230 uses the set of discrete material property parameters αi identified corresponding to the value of the discrete reference plastic strain ε ₀ ^p to identify the function-type parameters. A function indicating the relationship between the basic reference plastic strain ε ₀ ^p and the continuous material property parameter αi is identified (S34). In this case, a mathematical method such as a least square method is used.

  Using an aluminum alloy plate as a test material, the function of the material characteristic parameter αi of the yield function Yld2000-2d was identified. As a test material, an aluminum alloy plate 6016-T4 having a plate thickness of 0.9 mm was used. The main orientation is the Cube orientation, and the orientation density (random ratio) at a thickness of 1/4 is 133.

First, a uniaxial tensile test in the rolling direction is performed, and the uniaxial tensile stress σ _{0 at} the moment when the predetermined standard plastic strain ε ₀ ^p is reached and the plasticity per unit volume made until the standard plastic strain ε ₀ ^p is reached. It was determined the work _{W 0.} Next, a biaxial tensile test was performed, and a stress point (σ _x , σ _y ) when a plastic work equivalent to the plastic work W ₀ was consumed was obtained. They were plotted in the principal stress space to determine the isoplastic work surface relative to the reference plastic strain ε ₀ ^p .

FIG. 6 shows stresses σ _x and σ _y with respect to each standard plastic strain ε ₀ ^p . Although 13 standard plastic strains ε ₀ ^p were measured from 0.0002 to 0.16, 6 of them were shown as representatives in FIG.

  Next, using the simulation apparatus 200, the isoplastic work surface is approximated using the yield function Yld2000-2d. Here, Yld2000-2d is expressed by the following mathematical formulas 1 to 6, and has nine material characteristic parameters in total of order M and αi (i = 1 to 8).

Xi ′ and Xi ″ are the main values of the second floor tensors X ′ and X ″, and are calculated by the following equations 2 and 3.

X ′ and X ″ are obtained by linearly transforming the Cauchy stress tensor σ according to the following Equation 4.

The above Equation 4 can be displayed in matrix as Equations 5 and 6 below.

The parameter identification unit 230 has dimensionless uniaxial tensile stresses σ ₀ / σ ₀ , σ ₄₅ / σ ₀ , σ ₉₀ / σ ₀ and r values r ₀ , r in the 0 °, 45 °, and 90 ° directions from the rolling direction. _45, in addition to _{r 90,} determines the material properties parameter αi with equal biaxial tensile dimensionless stress was measured from the test σ _b / σ ₀ and plastic strain rate ratio ^{_{r b = dε y p / dε}} x p did. Further, the parameter identification unit 230 determines the order M so that the mean square error shown in the following Equation 7 is minimized with respect to the isoplastic work surface.

Here, n is the number of stress points constituting the non-dimensional plastic work surface, and k is the number of the linear stress path. r _k and r _{k 'is} the distance from the origin to yield surface that is actually measured stress points and modeled in the k-th linear stress path.

FIG. 7 shows a set of discrete material property parameter αi values respectively corresponding to the discrete reference plastic strain ε ₀ ^p values identified by the above calculation. It can be seen that the material characteristic parameter αi changes substantially continuously with respect to the change in the reference plastic strain ε ₀ ^p .

  The following formula 8 was selected as the function form of the order M.

The order M of Equation 8 was identified using the least square method. Using the order M identified in this way, the material property parameter αi corresponding to each reference plastic strain ε ₀ ^p was determined again.

The following formula 9 was selected as the function form of the material property parameter. Each parameter was identified by the least square method using the material characteristic parameter αi corresponding to each standard plastic strain ε ₀ ^p calculated above.

FIG. 8 shows parameters of each function identified as described above. In FIG. 7, the relationship between the material property parameter αi expressed by the function and the standard plastic strain ε ₀ ^p is shown by a curve. Furthermore, the result of approximating the isoplastic work surface using the material characteristic parameters M and αi represented by the function is shown by a curve in FIG. As can be seen from FIG. 6, it is understood that the work hardening anisotropy of the test material can be accurately reproduced by expressing the material characteristic parameters as M and αi as a function of the standard plastic strain ε ₀ ^p .

  FIG. 9 is a flowchart showing details of step S40 of FIG. 4 in the present embodiment. First, initialization of space discretization, yield function, etc. is performed (S41).

The deformation calculation unit 240 assigns a virtual time increment Δt (S42), and calculates an internal pressure increment ΔP corresponding to the time increment Δt (SS43). Further, the deformation calculation unit 240 gives a virtual displacement, that is, a strain increment to the object to be processed (S44). Deformation calculator 240, using a function representing the material properties parameter αi with reference plastic strain epsilon ₀ ^p, calculates the value of the material property parameters αi for the current value of the reference plastic strain ε ₀ ^p (S45). The deformation calculation unit 240 calculates the stress increment Δσ and the reference plastic strain increment Δε ₀ ^p using the value of the material property parameter αi (S46). In this case, the stress increment Δσ and the reference plastic strain increment Δε ₀ ^p are calculated based on the values of stress and strain in the repeated step n while maintaining material mechanical consistency. Based on the result of step S46, the deformation calculation unit 240 calculates a tangential stiffness coefficient (S47).

  The deformation calculation unit 240 confirms whether or not the virtual work principle can be satisfied with the virtual displacement given in step S44 (S48). If the condition cannot be satisfied, the process returns to step S44, and the calculation is repeated until the condition is satisfied (S48: No). When the condition is satisfied (S48: Yes), the deformation calculation unit 240 determines the virtual displacement and updates the coordinates (S49).

  Steps S42 to S49 are repeated with a preset number of repetition steps (S50). The deformation calculation unit 240 outputs the result of the deformation state when the repetition in step S50 is completed (S51).

In steps S45 and S46, every time the value of the material property parameter αi is used as the yield function, repeated calculation such as further correcting the value of the reference plastic strain ε ₀ ^p to satisfy the material mechanics condition is included. Good.

  An example of a specific method of steps S44 to S48 is as follows. First, the stress-strain relationship in which the plastic work increment is equivalent is defined as an equivalent stress-equivalent plastic strain diagram as shown in the following formula 10.

The backward Euler method is used for the stress integration algorithm. In order to express the anisotropic hardening behavior accompanying the increase in plastic work, an equivalent plastic strain ε- ^p is added to the argument. That is, in the yield function, the material characteristic parameter is introduced as a function of plastic strain. When the yield function is generally expressed, the following Expression 11 is obtained.

The direction of the plastic strain increment Δε ^p is expressed by the following Equation 12 when the related flow law is applied.

  Here, γ represents a plastic multiplier. The relationship between the plastic multiplier γ and the equivalent plastic strain increment is expressed by the following Equation 13 from Equations 11 and 12.

Here, since σ ⁻ (σ, ε− ^p ) is a linear function of σ, the relationship of the following formula 14 is used.

When the work hardening formula is expressed by the following formula 15, when the analysis step before the increment is n, the yield condition after the increment is expressed by the following formula 16, and the stress increment is expressed by the formula 17. Here, C ^e is the fourth floor of the elastic tensor.

  Equations 16 and 17 are solved using the Newton-Raphson method. In this case, multistage return mapping is used as a method for stably obtaining a converged solution. In this method, the yield function to be satisfied is relaxed by introducing substeps to maintain stable convergence, and the convergence is gradually removed to solve the convergence solution.

  In this case, the tangential stiffness is consistent with the stress integration algorithm. According to the backward Euler method, the yield condition and stress after increment, that is, n + 1 steps, satisfy the following formulas 18 and 19.

Taking σ _{n + 1} and ε _{n + 1} as variables and taking the total differentiation of both sides, the following formula 20 is obtained. Here, the definition of Equation 21 below was used.

As described above, anisotropic hardening of the specimen is introduced by expressing the material property parameter αi as a function of the reference plastic strain ε ₀ ^p in the simulation of the finite element method. As a result, the simulation reproduces the experiment more accurately.

  FIG. 10 schematically shows the hydraulic bulge forming used in Example 2. In order to verify the validity of the yield function, an experiment was verified by hydraulic bulge forming. A disk specimen 20 having a die hole diameter of φ150 mm was hydraulically bulged. The strain distribution was detected by a digital image correlation method using a pair of cameras 140 and 142.

  Abaqus / StandardVer 6.11-1 was used for the analysis. A mold having the same shape as the hydraulic bulge test was defined as an analytical rigid body, and the internal pressure was applied with the mold restraining force calculated by the bolt tightening force applied. The friction coefficient of the mold contact surface was 0.3. In consideration of symmetry, a 1/4 model was used, and a 4-node reduced integral shell element was used. Element division is 1 mm in the radial direction (however, 0.5 mm in the range of bead-formed radius 92 to 103 mm), the range of radius 30 to 110 mm is 36 in the circular direction, and the range within 20 mm radius is about 1 mm per side Automatic division.

As the work hardening formula, the Swift formula identified by the rolling direction uniaxial tensile test was used. Yld2000-2d identified in Example 1 was used as the yield function used for the analysis. That is, Yld2000-2d was used in which each material characteristic parameter M, αi was expressed as a function of the plastic strain ε ^p identified by each parameter in Equations 8 and 9 and FIG. Hereinafter, this yield function may be expressed as Yld2000-2d (DWH) or simply DWH. The finite element method was executed by the deformation calculation unit 240 using the above conditions.

As comparative examples, von Mises, which is an isotropic hardening model, and Hill'48 yield quadratic function were used. As another comparative example, Yld2000-2d was used in which the material property parameters M and αi were identified with respect to the isoplastic work surface with respect to the standard plastic strain ε ₀ ^p of 0.002, 0.04, and 0.16. These material characteristic parameters αi are shown in FIG.

FIG. 12 shows a comparison between the experimental value of the true plate thickness strain | ε _z | The end point of the finite element method calculation indicates that the calculation does not converge unless the internal pressure increment is set to ΔP = 1 × 10 ⁻¹⁰ MPa / step or less, and corresponds to the maximum internal pressure applied in the experiment.

  As shown in FIG. 12, Yld2000-2d (DWH) substantially matches the experimental value over the entire region, and the prediction accuracy of the maximum internal pressure is also improved. In contrast, von Mises and M = 6.7 overestimate the maximum internal pressure. On the other hand, Hill'48 predicts the maximum internal pressure too low. For M = 12.5 and 37.8, the experimental value is well reproduced when the internal pressure is 4 MPa or less, but the maximum internal pressure is predicted excessively.

  FIGS. 13 to 15 are thickness strain distributions along meridians in the 0 °, 45 °, and 90 ° directions from the rolling direction. FIG. 13 shows the bulge forming height h = 10 mm, FIG. 14 shows 20 mm, and FIG. 15 shows 38 mm.

In FIG. 13, the plate thickness strain | ε _z | at the apex is about 0.025, which corresponds to the low strain state. Yld2000-2d (DWH) can reproduce the experimental values most accurately over the entire area. On the other hand, although M = 6.7 can reproduce the tendency of the experimental value well, other vonMises, Hill'48, M = 12.5, 38.7 are overvalued in any region. The result is underestimated.

In FIG. 14, the plate thickness strain | ε _z | at the apex is about 0.08, which corresponds to a moderate strain state. Yld2000-2d (DWH) can reproduce the experimental values most accurately over the entire area. On the other hand, although M = 12.5 can reproduce the tendency of the experimental value well, other vonMises, Hill'48, M = 6.7, 38.7 are overvalued in any region. The result is underestimated.

In FIG. 15, the plate thickness strain | ε _z | at the apex is about 0.28, which corresponds to a large strain state. Yld2000-2d (DWH) can accurately reproduce the experimental values in the range of R from 0 to 20 mm. On the other hand, von Mises, Hill'48, M = 6.7, 12.5, and 38.7 are the results of evaluating the value in either region as being too large or too small.

  As described above, according to the results of FIGS. 12 to 15, all the moldings from the initial stage of molding to the end of molding are performed by using Yld2000-2d in which the material characteristic parameters are expressed as a function of plastic strain as Example 2. Analysis accuracy improved for each stage.

  As described above, according to the present embodiment, by expressing the material characteristic parameter as a function of the standard plastic strain, it is possible to perform a simulation that more accurately reproduces anisotropic hardening. Until now, even though actual processing often involves anisotropic hardening, the reproducibility of anisotropic hardening has been insufficient in simulations, so large-scale destructive testing of the actual product has been performed. On the other hand, according to the present embodiment, the deformation state of the object to be processed in consideration of anisotropic hardening can be predicted, so that a destructive test or the like is not required, and time and cost can be greatly reduced. Moreover, since the information used when identifying the material characteristic parameter with respect to the standard plastic strain is obtained from the biaxial bulge test apparatus having a variable stress ratio, the more accurate material parameter can be specified for each material.

  As mentioned above, although this invention was demonstrated using embodiment, the technical scope of this invention is not limited to the range as described in the said embodiment. It will be apparent to those skilled in the art that various modifications or improvements can be added to the above-described embodiment. It is apparent from the scope of the claims that the embodiments added with such changes or improvements can be included in the technical scope of the present invention.

  The order of execution of each process such as operations, procedures, steps, and stages in the apparatus, system, program, and method shown in the claims, the description, and the drawings is particularly “before” or “prior to”. It should be noted that the output can be realized in any order unless the output of the previous process is used in the subsequent process. Regarding the operation flow in the claims, the description, and the drawings, even if it is described using “first”, “next”, etc. for convenience, it means that it is essential to carry out in this order. It is not a thing.

  10 disc test piece, 12 hollow part, 20 disc test piece, 100 biaxial bulge test device, 102 body, 104 controller, 110, 112 gripping part, 120, 122 driving part, 130 pressure adjusting part, 132 flow path, 140 , 142 camera, 200 simulation device, 202 CD-ROM, 210 parameter acquisition unit, 220 stress value acquisition unit, 230 parameter identification unit, 240 deformation calculation unit, 250 storage unit

Claims (7)

A simulation device that performs a molding simulation on an object of a predetermined material using a predetermined yield function including a material property parameter,
For the material, a parameter acquisition unit for acquiring a relationship between a value of a standard plastic strain and a value of the material property parameter corresponding to the standard plastic strain;
When given a virtual displacement to said object, successive values of the material properties parameter corresponding to the amount of standard plastic strain according to the virtual time increment for each incremental the virtual time, using the yield function And a deformation calculation unit that calculates a deformation state of the object. The parameter acquisition unit
A plurality of biaxial stress values when a biaxial tensile test is performed on the material so that the biaxial stress ratios are different are obtained in association with a plurality of predetermined reference plastic strain values. A stress value acquisition unit;
The simulation apparatus according to claim 1, further comprising: a parameter identification unit that identifies a value of the material property parameter for each of the plurality of reference plastic strain values based on the plurality of biaxial stress values.   The parameter identification unit identifies the value of the material characteristic parameter so that a yield surface for each of the plurality of reference plastic strain values matches the stress value of the isoplastic work surface in the yield function. The simulation apparatus described in 1.   The said parameter identification part identifies the function which shows the relationship between a reference | standard plastic strain and a material characteristic parameter based on the value of the said material characteristic parameter identified for every value of these some reference | standard plastic strain. Simulation device.   The simulation apparatus according to claim 3, wherein the deformation calculation unit further corrects the value of the reference plastic strain so as to satisfy a condition of material mechanics every time the value of the material characteristic parameter is used for the yield function. A simulation method for performing a molding simulation on an object of a predetermined material using a predetermined yield function including a material property parameter,
For the material, a parameter acquisition step of acquiring a relationship between a value of a standard plastic strain and a value of the material property parameter corresponding to the standard plastic strain;
When given a virtual displacement to said object, successive values of the material properties parameter corresponding to the amount of standard plastic strain according to the virtual time increment for each incremental the virtual time, using the yield function A deformation calculation step of calculating a deformation state of the object. A program for controlling a computer that performs a molding simulation on an object of a predetermined material using a predetermined yield function including a material property parameter,
For the material, a parameter acquisition step of acquiring a relationship between a value of a standard plastic strain and a value of the material property parameter corresponding to the standard plastic strain;
When given a virtual displacement to said object, successive values of the material properties parameter corresponding to the amount of standard plastic strain according to the virtual time increment for each incremental the virtual time, using the yield function And a deformation calculation step of calculating a deformation state of the object.