{"id":354,"date":"2022-07-29T13:19:40","date_gmt":"2022-07-29T13:19:40","guid":{"rendered":"http:\/\/summergeometry.org\/sgi2022\/?p=354"},"modified":"2022-07-29T13:26:40","modified_gmt":"2022-07-29T13:26:40","slug":"using-gray-scott-reaction-diffusion-model-to-create-spots-and-stripes-on-3d-spot","status":"publish","type":"post","link":"https:\/\/summergeometry.org\/sgi2022\/using-gray-scott-reaction-diffusion-model-to-create-spots-and-stripes-on-3d-spot\/","title":{"rendered":"Using Gray-Scott Reaction-Diffusion Model to Create Spots (and Stripes) on 3D Spot"},"content":{"rendered":"\n

By SGI Fellow Natalie Patten<\/em><\/p>\n\n\n\n

The Gray-Scott reaction-diffusion model consists of two partial differential equations: $$\\frac{\\partial U}{\\partial t}=D_u\\nabla^2U-UV^2+F(1-U)$$ $$\\frac{\\partial V}{\\partial t}=D_v\\nabla^2V+UV^2-(F+k)V.$$Using a semi-implicit approach, each \\(U^{i+1}\\) can be calculated as follows: \\((I -{\\Delta}tD_uL)U^{i+1}=U^i + {\\Delta}t (F(1-U^i) – U^iV^iV^i)\\). Here, \\(I\\) is the identity matrix, \\({\\Delta}t\\) is change in time, \\(D_u\\) is the diffusion rate, \\(L\\) is the Laplacian of the mesh, \\(F\\) is the feed rate, and \\(k\\) is the degrading rate. Similarly, \\(V^{i+1}\\) can be calculated as follows: \\((I -{\\Delta}tD_vL)V^{i+1}=V^i + {\\Delta}t (U^iV^iV^i-(F+k)V^i)\\). Using Eigen<\/a>‘s SimplicialLDLT solver, which uses Cholesky decomposition, we can solve for each: <\/strong><\/p>\n\n\n

\nEigen::SparseMatrix<double> temp1 = (I - timeStep * du * L).eval();\nEigen::SparseMatrix<double> temp2 = (I - timeStep * dv * L).eval();\n\nEigen::SimplicialLDLT<Eigen::SparseMatrix<double>> solver1;\nEigen::SimplicialLDLT<Eigen::SparseMatrix<double>> solver2;\n\nsolver1.compute(temp1);\nsolver2.compute(temp2);\n\nfor (int i = 0; i < numSteps; i++) {\n   Eigen::VectorXd UVV = U.array() * V.array() * V.array();\n   Eigen::VectorXd FU = F.array() * U.array();\n   U = solver1.solve((U + timeStep * (F - FU - UVV)).eval());\n   Eigen::VectorXd kV = k.array() * V.array();\n   Eigen::VectorXd FV = F.array() * V.array();\n   V = solver2.solve((V + timeStep * (UVV - kV - FV)).eval());\n}\n<\/pre><\/div>\n\n\n
U was initialized to be 1.0 everywhere, and V was initialized to be 0.0 everywhere except in one small area where it was initialized to be 1.0. Du<\/sub> and Dv<\/sub> were set to be 1.0 and 0.5, respectively. I then modelled this onto a 3D mesh of Spot the cow using Polyscope<\/a> with a time step of 1 and various F and k values, based on parameters found in Pearson’s “Complex Patterns in a Simple System” (1993). Below is a summary of the initial results. You can see various patterns formed by the Gray-Scott model.<\/p>\n\n\n\n
F = 0.024 and k = 0.06 at 10,000 steps:<\/h2>\n\n\n\n
\"\"<\/figure>\n\n\n\n
F = 0.04 and k = 0.06 at 10,000 steps:<\/h2>\n\n\n\n
\"\"<\/figure>\n\n\n\n
F = 0.05 and k = 0.06 at 10,000 steps:<\/h2>\n\n\n\n
\"\"<\/figure>\n\n\n\n
F = 0.03 and k = 0.056 at 10,000 steps:<\/h2>\n\n\n\n
\"\"<\/figure>\n\n\n\n
F = 0.044 and k =  0.063 at 10,000 steps:<\/h2>\n\n\n\n
\"\"<\/figure>\n\n\n\n
Acknowledgments:<\/p>\n\n\n\n
  • Thank you to my project partners Mariem Khlifi and Gimin Nam for help with (lots of!) debugging.<\/li>
  • Thank you to my project mentor Etienne Vouga and project TA Erick Jimenez Berumen for guidance, debugging, and resources.<\/li><\/ul>\n","protected":false},"excerpt":{"rendered":"
    By SGI Fellow Natalie Patten The Gray-Scott reaction-diffusion model consists of two partial differential equations: $$\\frac{\\partial U}{\\partial t}=D_u\\nabla^2U-UV^2+F(1-U)$$ $$\\frac{\\partial V}{\\partial t}=D_v\\nabla^2V+UV^2-(F+k)V.$$Using a semi-implicit approach, each \\(U^{i+1}\\) can be calculated as follows: \\((I -{\\Delta}tD_uL)U^{i+1}=U^i + {\\Delta}t (F(1-U^i) – U^iV^iV^i)\\). Here, \\(I\\) is the identity matrix, \\({\\Delta}t\\) is change in time, \\(D_u\\) is the diffusion rate, \\(L\\) […]<\/p>\n","protected":false},"author":29,"featured_media":0,"comment_status":"open","ping_status":"open","sticky":false,"template":"","format":"standard","meta":{"footnotes":""},"categories":[45],"tags":[40],"class_list":["post-354","post","type-post","status-publish","format-standard","hentry","category-research","tag-reaction-diffusion"],"_links":{"self":[{"href":"https:\/\/summergeometry.org\/sgi2022\/wp-json\/wp\/v2\/posts\/354","targetHints":{"allow":["GET"]}}],"collection":[{"href":"https:\/\/summergeometry.org\/sgi2022\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/summergeometry.org\/sgi2022\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/summergeometry.org\/sgi2022\/wp-json\/wp\/v2\/users\/29"}],"replies":[{"embeddable":true,"href":"https:\/\/summergeometry.org\/sgi2022\/wp-json\/wp\/v2\/comments?post=354"}],"version-history":[{"count":9,"href":"https:\/\/summergeometry.org\/sgi2022\/wp-json\/wp\/v2\/posts\/354\/revisions"}],"predecessor-version":[{"id":556,"href":"https:\/\/summergeometry.org\/sgi2022\/wp-json\/wp\/v2\/posts\/354\/revisions\/556"}],"wp:attachment":[{"href":"https:\/\/summergeometry.org\/sgi2022\/wp-json\/wp\/v2\/media?parent=354"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/summergeometry.org\/sgi2022\/wp-json\/wp\/v2\/categories?post=354"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/summergeometry.org\/sgi2022\/wp-json\/wp\/v2\/tags?post=354"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}