[Insight-users] Implementation Finite Element Method (FEM) in ITK

[Jeremy Bournesel]

Thanks Brian. you comments are helping me tremendously.

> we use standard gaussian integration points to integrate the
> derivatives of an image similarity over each element and distribute
> these derivatives to the nodes of the elements.

Where is this actually happening in the code. Is there any kind of
documentation about this.

> not sure exactly what you mean here, but the elements are placed into
> the global linear system that models the whole problem in the physical
> space.

When and how is the global linear system created?

Again, I really appreciate your input.

-J


On Wed, Feb 2, 2011 at 10:19 PM, brian avants <stnava at gmail.com> wrote:

> jeremy
>
> sorry, i was really responding to a question cristina asked some time
> ago and it took me a bit to remember the relevant pieces of math to
> properly answer her .
>
> > Are you able to explain the bigger picture?
>
> in fem registration, you want to "deform" an image slowly over time
> until the external image forces reach equilibrium with the internal
> body forces (the regularization provided by the physical model).  the
> parameters of the model control all of the above interactions.
>
> > I have the reference and the target image.
> > If no grid/mesh is specified a uniform mesh is created.
>
> yes, i think that's correct.
>
> > I assume there are no initial forces. But what happens then?
>
> there is one fixed boundary condition to make the problem well-posed.
> otherwise, no forces.
>
> > How are the forces and deformations for every node in every element
> > calculated?
>
> we use standard gaussian integration points to integrate the
> derivatives of an image similarity over each element and distribute
> these derivatives to the nodes of the elements.
>
> > Is it based on pixel data (histogram matching) close to the nodes?
>
> there are a few options but you've got the basic idea --- you're
> trying to match the local intensity.
>
> > Based on the local energy at the element level, how is the global energy
> > calculated?
>
> not sure exactly what you mean here, but the elements are placed into
> the global linear system that models the whole problem in the physical
> space.
>
> b.a.
>
> > These are the kind of questions I have, it be really awesome if you can
> > help.
> > Thanks,
> > -J
> >
> >
> > On Tue, Feb 1, 2011 at 8:17 PM, brian avants <stnava at gmail.com> wrote:
> >>
> >> hi jeremy , cristina
> >>
> >> here is a brief latex sketch of the method used in the
> >> FEMCrankNicolsonSolver --
> >>
> >> \text{matrix form parabolic PDE, heat equation type} \\
> >> M\frac{ \partial U}{\partial t} +  K U  = f \\
> >> \frac{ M U_t - M  U_{t-1}}{\delta_t} +  K U  = f \\
> >> \frac{ M U_t - M  U_{t-1} +  \delta_t K U }{\delta_t} = f \\
> >> \text{ here i make the decision to represent $U$ and $f$ at $t$ rather
> >> than $t-1$ i believe this is backward Euler}\\
> >> ( M  + \delta_t K ) U_t  =  M  U_{t-1} + \delta_t f_t \\
> >> \\
> >> \text{other choices are possible, e.g. } \\
> >> \text{ Crank-Nicolson discretization is based on averaging forward and
> >> backward Euler } \\
> >> ( M + \alpha \delta_t K ) U_t = (  M - ( 1- \alpha ) \delta_t  K )
> >> U_{t-1} + f \\
> >>  f  =  \delta_t ( \alpha f_t + ( 1 - \alpha ) f_{t-1} ) \\
> >> \\
> >> \text{if}~ \alpha = 1 \text{ backward Euler }\\
> >> ( M +  \delta_t K ) U_t = M  U_{t-1} + \delta_t f \\
> >> \\
> >> \text{if}~ \alpha = 0 \text{ forward Euler }\\
> >>  M U_t =(  M - \delta_t  K ) U_{t-1} + \delta_t f_{t-1}
> >>
> >> hopefully this is helpful.
> >>
> >> brian
> >>
> >>
> >>
> >>
> >> On Tue, Feb 1, 2011 at 4:01 PM, Jeremy Bournesel
> >> <jeremy.bournesel at gmail.com> wrote:
> >> > Hi,
> >> > I have a decent background regarding FEM in the engineering domain and
> >> > have
> >> > started looking at it for registering medical images as well.
> >> > I ran one of the ITK deformable registration demos, read the Software
> >> > Guide
> >> > and the ITK Powerpoint Presentations, however there is still some
> magic
> >> > going on that I don't completely get.
> >> > From what I understand ITK is using a uniform grid/mesh (if a custom
> one
> >> > wasn't supplied) and then calculates the deformation for every element
> >> > at
> >> > the nodal points (iteratively).
> >> > I couldn't find any information on how the exact process is working
> >> > (besides
> >> > some high level slides of the type K U = F).
> >> > It'd be awesome if someone can point me to the right document or
> explain
> >> > it
> >> > to me.
> >> > Thanks,
> >> > Jeremy
> >> >
> >> >
> >> >
> >> >
> >> >
> >> > _____________________________________
> >> > Powered by www.kitware.com
> >> >
> >> > Visit other Kitware open-source projects at
> >> > http://www.kitware.com/opensource/opensource.html
> >> >
> >> > Kitware offers ITK Training Courses, for more information visit:
> >> > http://www.kitware.com/products/protraining.html
> >> >
> >> > Please keep messages on-topic and check the ITK FAQ at:
> >> > http://www.itk.org/Wiki/ITK_FAQ
> >> >
> >> > Follow this link to subscribe/unsubscribe:
> >> > http://www.itk.org/mailman/listinfo/insight-users
> >> >
> >> >
> >
> >
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20110203/2b48b216/attachment.htm>