Systems, methods and, media for simulating deformations of nonlinear elastic bodies

1. A method for simulating deformation of an elastic body, the method comprising:

determining, using a hardware processor, for each of a plurality of macroblocks B i , including a first macroblock B 1 , a stiffness matrix K i corresponding to at least a portion of a model of a non-linear elastic solid that is partitioned into a plurality of cells, wherein entries in the stiffness matrix K i correspond to nodes associated with cells of the macroblock B i ;

converting, for the first macroblock B 1 , the stiffness matrix K 1 into block form to include a submatrix K I 1 I 1 , where the subscript I 1 I 1 indicates that the submatrix K I 1 I 1 includes entries in which both nodes associated with the entry are between cells of the first macroblock B 1 ;

determining, for the first macroblock B 1 , at least a portion of an inverse matrix K I 1 I 1 −1 of the submatrix K I 1 I 1 ;

receiving input data corresponding to force applied to one or more nodes of the plurality of cells of the model;

determining, for the first macroblock B 1 , displacements of nodes that are on the exterior of the first macroblock B 1 based at least in part on the input data and the portion of the inverse matrix K I 1 I 1 −1 ;

determining, for the first macroblock B 1 , displacements of nodes that are interior to the first macroblock B 1 based at least in part on the input data and the displacements of nodes that are on the exterior of the first macroblock B 1 ;

determining updated positions of the cells of the model based at least in part on the displacements of nodes that are on the exterior of the macroblocks; and

causing the model to be presented on a display device using the updated positions.

2. The method of claim 1 , wherein each of the plurality of macroblocks B i has 16×8×8 grid cells and comprise 15×7×7 internal nodes, wherein each macroblock is associated with a submatrix K I i I i that only includes entries for nodes between two cells of macroblock B i .

3. The method of claim 2 , wherein determining at least the portion of the inverse matrix K I 1 I 1 −1 further comprises:

partitioning the 15×7×7 internal nodes of the first macroblock B 1 into sixteen 3×3×3 subdomains and five interface layers, wherein a first interface layer represents a 1×7×7 layer of nodes that separates the 15×7×7 internal nodes into a first 7×7×7 subdomain and a second 7×7×7 subdomain;

generating a block form of the submatrix K I 1 I 1 by reordering the submatrix K I 1 I 1 into a first submatrix K 11 corresponding to entries representing nodes in the first 7×7×7 subdomain, a second submatrix K 22 corresponding to entries representing nodes in the second 7×7×7 subdomain, a third submatrix K cc corresponding to entries representing nodes in the first interface layer, a submatrix K 1c corresponding to entries representing nodes in the first 7×7×7 subdomain and the first interface layer, a submatrix K 2c corresponding to entries between nodes in the second 7×7×7 subdomain and the first interface layer, a submatrix K c1 corresponding to entries between nodes in the first interface layer and the first 7×7×7 subdomain, and a submatrix K c2 corresponding to entries between nodes in the first interface layer and the second 7×7×7 subdomain;

generating at least a portion of the inverse of the block form of the submatrix K I 1 I 1 by converting the block form of the submatrix K I 1 I 1 to a block-LDL form:

(

I

0

-

K

11

-

1

⁢

K

1

⁢

c

0

I

-

K

22

-

1

⁢

K

2

⁢

c

0

0

I

)

⁢

(

K

11

-

1

0

0

0

K

22

-

1

0

0

0

C

)

⁢

(

I

0

0

0

I

0

-

K

c

⁢

⁢

1

⁢

K

11

-

1

-

K

c

⁢

⁢

2

⁢

K

22

-

1

I

)

where C in the block-LDL form represents the Shur complement of K cc and is equal to K cc −K c1 K 11 −1 K 1c −K c2 K 22 −1 K 2c ;

determining, for each of the sixteen 3×3×3 subdomains, an inverse K jj −1 of a submatrix K jj corresponding to entries within that 3×3×3 subdomain, where subscript j represents which of the sixteen 3×3×3 subdomains the submatrix represents;

determining K 11 −1 based on each inverse K jj −1 that corresponds to a 3×3×3 subdomain included within the first subdomain;

determining K 22 −1 based on each inverse K jj −1 that corresponds to a 3×3×3 subdomain included within the second subdomain; and

storing at least the portion of the inverse matrix K I 1 I 1 −1 in cache memory of the hardware processor.

4. The method of claim 3 , wherein inverting each submatrix K jj further comprises:

determining Cholesky factors for the submatrix K jj ; and

determining the inverse K jj −1 using the Cholesky factors and forward and backward substitution.

5. The method of claim 3 , further comprising inverting the submatrices K jj corresponding to the sixteen 3×3×3 subdomains in parallel.

6. The method of claim 1 , further comprising determining at least the portion of an inverse matrix K I i I i −1 corresponding to each of the plurality of macroblocks B i , including the portion of the inverse matrix K I 1 I 1 −1 corresponding the first macroblock B 1 , in parallel.

7. The method of claim 1 , wherein determining displacements of nodes that are on the exterior of the first macroblock B 1 based at least in part on the input data and at least the portion of the inverse matrix K I 1 I 1 −1 further comprises using a Conjugate Gradient-based iterative solver.

8. The method of claim 1 , further comprising:

determining, for each of a second plurality of macroblocks B i that each overlap a boundary of the model, a stiffness matrix K i corresponding to at least a portion of the model, wherein each macroblock of the second plurality of macroblocks includes at least one empty cell; and

setting, for each of the second plurality of macroblocks, entries in the stiffness matrix K i that correspond to nodes in the macroblock that are exterior to the model to a zero-Dirichlet condition.

9. A system for simulating deformation of an elastic body, the system comprising:

memory storing a model of a non-linear elastic solid;

a display device;

a hardware processor that is coupled to the memory and the display device, and is programmed to;

determine, using a hardware processor, for each of a plurality of macroblocks B i , including a first macroblock B 1 , a stiffness matrix K i corresponding to at least a portion of the model that is partitioned into a plurality of cells, wherein entries in the stiffness matrix correspond to nodes associated with cells of the macroblock B i ;

convert, for for first macroblock B 1 , the stiffness matrix K 1 into block form to include a submatrix K I 1 I 1 , where the subscript I 1 I 1 indicates that the submatrix K I 1 I 1 includes entries in which both nodes associated with the entry are between cells of the first macroblock B 1 ;

determining, for the first macroblock B 1 , at least a portion of an inverse matrix K I 1 I 1 −1 of the submatrix K I 1 I 1 ;

receive input data corresponding to force applied to one or more nodes of the plurality of cells of the model;

determine, for the first macroblock B 1 , displacements of nodes that are on the exterior of the first macroblock B 1 based at least in part on the input data and the portion of the inverse matrix K I 1 I 1 −1 ;

determine, for the first macroblock B 1 , displacements of nodes that are interior to the first macroblock B 1 based at least in part on the input data and the displacements of nodes that are on the exterior of the first macroblock B 1 ;

determine updated positions of the cells of the model based at least in part on the displacements of nodes that are on the exterior of the macroblocks; and

cause the model to be presented on the display device using the updated positions.

10. The system of claim 9 , wherein each of the plurality of macroblocks B i has 16×8×8 grid cells and comprise 15×7×7 internal nodes, wherein each macroblock is associated with a submatrix K I i I i that only includes entries for nodes between two cells of macroblock B i .

11. The system of claim 10 , wherein the hardware processor is associated with cache memory, and is further programmed to:

partition the 15×7×7 internal nodes of the first macroblock B 1 into sixteen 3×3×3 subdomains and five interface layers, wherein a first interface layer represents a 1×7×7 layer of nodes that separates the 15×7×7 internal nodes into a first 7×7×7 subdomain and a second 7×7×7 subdomain;

generate a block form of the submatrix K I 1 I 1 by reordering the submatrix K I 1 I 1 into a first submatrix K 11 corresponding to entries representing nodes in the first 7×7×7 subdomain, a second submatrix K 22 corresponding to entries representing nodes in the second 7×7×7 subdomain, a third submatrix K cc corresponding to entries representing nodes in the first interface layer, a submatrix K 1c corresponding to entries representing nodes in the first 7×7×7 subdomain and the first interface layer, a submatrix K 2c corresponding to entries between nodes in the second 7×7×7 subdomain and the first interface layer, a submatrix K c1 corresponding to entries between nodes in the first interface layer and the first 7×7×7 subdomain, and a submatrix K c2 corresponding to entries between nodes in the first interface layer and the second 7×7×7 subdomain;

generate at least a portion of the inverse of the block form of the submatrix K I 1 I 1 by converting the block form of the submatrix K I 1 I 1 to a block-LDL form:

(

I

0

-

K

11

-

1

⁢

K

1

⁢

c

0

I

-

K

22

-

1

⁢

K

2

⁢

c

0

0

I

)

⁢

(

K

11

-

1

0

0

0

K

22

-

1

0

0

0

C

)

⁢

(

I

0

0

0

I

0

-

K

c

⁢

⁢

1

⁢

K

11

-

1

-

K

c

⁢

⁢

2

⁢

K

22

-

1

I

)

where C in the block-LDL form represents the Shur complement of K cc and is equal to K cc −K c1 K 11 −1 K 1c −K c2 K 22 −1 K 2c ;

determine, for each of the sixteen 3×3×3 subdomains, an inverse K jj −1 of a submatrix K jj corresponding to entries within that 3×3×3 subdomain, where subscript j represents which of the sixteen 3×3×3 subdomains the submatrix represents;

determine K 11 −1 based on each inverse K jj −1 that corresponds to a 3×3×3 subdomain included within the first subdomain;

determine K 22 −1 based on each inverse K jj −1 that corresponds to a 3×3×3 subdomain included within the second subdomain; and

store at least the portion of the inverse matrix K I 1 I 1 −1 in cache memory of the hardware processor.

12. The system of claim 11 , wherein the hardware processor is further programmed to:

determine, for each submatrix K jj , Cholesky factors for the submatrix K jj ; and

determine, for each submatrix K jj , the inverse K jj −1 using the Cholesky factors and forward and backward substitution.

13. The system of claim 11 , wherein the hardware processor is further programmed to invert the submatrices K jj corresponding to the sixteen 3×3×3 subdomains in parallel using single instruction, multiple data instructions.

14. The system of claim 9 , wherein the hardware processor is further programmed to determine at least the portion of an inverse matrix K I i I i −1 corresponding to each of the plurality of macroblocks B i , including the portion of the inverse matrix K I 1 I 1 −1 corresponding the first macroblock B 1 , in parallel using single instruction, multiple data instructions.

15. The system of claim 9 , wherein the hardware processor is further configured to use a Conjugate Gradient-based iterative solver to determine displacements of nodes that are on the exterior of the first macroblock B 1 .

16. The system of claim 9 , wherein the hardware processor is further programmed to:

determine, for each of a second plurality of macroblocks B i that each overlap a boundary of the model, a stiffness matrix K i corresponding to at least a portion of the model, wherein each macroblock of the second plurality of macroblocks includes at least one empty cell; and

set, for each of the second plurality of macroblocks, entries in the stiffness matrix K i that correspond to nodes in the macroblock that are exterior to the model to a zero-Dirichlet condition.

17. A non-transitory computer readable medium containing computer executable instructions that, when executed by a processor, cause the processor to perform a method for simulating deformation of an elastic body, the method comprising:

determining for each of a plurality of macroblocks B i , including a first macroblock B 1 , a stiffness matrix K i corresponding to at least a portion of a model of a non-linear elastic solid that is partitioned into a plurality of cells, wherein entries in the stiffness matrix K i correspond to nodes associated with cells of the macroblock B i ;

converting, for the first macroblock B 1 , the stiffness matrix K 1 into block form to include a submatrix K I 1 I 1 , where the subscript I 1 I 1 indicates that the submatrix K I 1 I 1 includes entries in which both nodes associated with the entry are between cells of the first macroblock B 1 ;

determining, for the first macroblock B 1 , at least a portion of an inverse matrix K I 1 I 1 −1 of the submatrix K I 1 I 1 ;

receiving input data corresponding to force applied to one or more nodes of the plurality of cells of the model;

determining, for the first macroblock B 1 , displacements of nodes that are on the exterior of the first macroblock B 1 based at least in part on the input data and the portion of the inverse matrix K I 1 I 1 −1 ;

determining, for the first macroblock B 1 , displacements of nodes that are interior to the first macroblock based B 1 at least in part on the input data and the displacements of nodes that are on the exterior of the first macroblock B 1 ;

determining updated positions of the cells of the model based at least in part on the displacements of nodes that are on the exterior of the macroblocks; and

causing the model to be presented on a display device using the updated positions.

18. The non-transitory computer readable medium of claim 17 , wherein each of the plurality of macroblocks B i has 16×8×8 grid cells and comprise 15×7×7 internal nodes, wherein each macroblock is associated with a submatrix K I i I i that only includes entries for nodes between two cells of macroblock B i .

19. The non-transitory computer readable medium of claim 18 , wherein determining at least the portion of the inverse matrix K I 1 I 1 −1 further comprises:

partitioning the 15×7×7 internal nodes of the first macroblock B 1 into sixteen 3×3×3 subdomains and five interface layers, wherein a first interface layer represents a 1×7×7 layer of nodes is that separates the 15×7×7 internal nodes into a first 7×7×7 subdomain and a second 7×7×7 sub domain;

generating a block form of the submatrix K I 1 I 1 by reordering the submatrix K I 1 I 1 into a first submatrix K 11 corresponding to entries representing nodes in the first 7×7×7 subdomain, a second submatrix K 22 corresponding to entries representing nodes in the second 7×7×7 subdomain, a third submatrix K cc corresponding to entries representing nodes in the first interface layer, a submatrix K 1c corresponding to entries representing nodes in the first 7×7×7 subdomain and the first interface layer, a submatrix K 2c corresponding to entries between nodes in the second 7×7×7 subdomain and the first interface layer, a submatrix K c1 corresponding to entries between nodes in the first interface layer and the first 7×7×7 subdomain, and a submatrix K c2 corresponding to entries between nodes in the first interface layer and the second 7×7×7 subdomain;

generating at least a portion of the inverse of the block form of the submatrix K I 1 I 1 by converting the block form of the submatrix K I 1 I 1 to a block-LDL form:

(

I

0

-

K

11

-

1

⁢

K

1

⁢

c

0

I

-

K

22

-

1

⁢

K

2

⁢

c

0

0

I

)

⁢

(

K

11

-

1

0

0

0

K

22

-

1

0

0

0

C

)

⁢

(

I

0

0

0

I

0

-

K

c

⁢

⁢

1

⁢

K

11

-

1

-

K

c

⁢

⁢

2

⁢

K

22

-

1

I

)

where C in the block-LDL form represents the Shur complement of K cc and is equal to K cc −K c1 K 11 −1 K 1c −K c2 K 22 −1 K 2c ;

determining, for each of the sixteen 3×3×3 subdomains, an inverse K jj −1 of a submatrix K jj corresponding to entries within that 3×3×3 subdomain, where subscript j represents which of the sixteen 3×3×3 subdomains the submatrix represents;

determining K 11 −1 based on each inverse K jj −1 that corresponds to a 3×3×3 subdomain included within the first subdomain;

determining K 22 −1 based on each inverse K jj −1 that corresponds to a 3×3×3 subdomain included within the second subdomain; and

storing at least the portion of the inverse matrix K I 1 I 1 −1 in cache memory of the hardware processor.

20. The non-transitory computer readable medium of claim 19 , wherein inverting each submatrix K jj further comprises:

determining Cholesky factors for the submatrix K jj ; and

determining the inverse K jj −1 using the Cholesky factors and forward and backward substitution.

21. The non-transitory computer readable medium of claim 19 , wherein the method further comprises inverting the submatrices K jj corresponding to the sixteen 3×3×3 subdomains in parallel.

22. The non-transitory computer readable medium of claim 17 , wherein the method further comprises determining at least the portion of an inverse matrix K I i I i −1 corresponding to each of the plurality of macroblocks B i , including the portion of the inverse matrix K I 1 I 1 −1 corresponding the first macroblock B 1 , in parallel.

23. The non-transitory computer readable medium of claim 17 , wherein determining displacements of nodes that are on the exterior of the first macroblock B 1 based at least in part on the input data and at least the portion of the inverse matrix K I 1 I 1 −1 further comprises using a Conjugate Gradient-based iterative solver.

24. The non-transitory computer readable medium of claim 17 , wherein the method further comprises:

determining, for each of a second plurality of macroblocks B i that each overlap a boundary of the model, a stiffness matrix K i corresponding to at least a portion of the model, wherein each macroblock of the second plurality of macroblocks includes at least one empty cell; and

setting, for each of the second plurality of macroblocks, entries in the stiffness matrix K i that correspond to nodes in the macroblock that are exterior to the model to a zero-Dirichlet condition.