Rocstar  1.0
Rocstar multiphysics simulation application
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
determinant.h File Reference
Include dependency graph for determinant.h:
This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Functions

template<class FT >
CGAL_BEGIN_NAMESPACE FT det2x2_by_formula (const FT &a00, const FT &a01, const FT &a10, const FT &a11)
 
template<class FT >
CGAL_KERNEL_MEDIUM_INLINE FT det3x3_by_formula (const FT &a00, const FT &a01, const FT &a02, const FT &a10, const FT &a11, const FT &a12, const FT &a20, const FT &a21, const FT &a22)
 
template<class FT >
CGAL_KERNEL_LARGE_INLINE FT det4x4_by_formula (const FT &a00, const FT &a01, const FT &a02, const FT &a03, const FT &a10, const FT &a11, const FT &a12, const FT &a13, const FT &a20, const FT &a21, const FT &a22, const FT &a23, const FT &a30, const FT &a31, const FT &a32, const FT &a33)
 
template<class FT >
CGAL_KERNEL_LARGE_INLINE FT det5x5_by_formula (const FT &a00, const FT &a01, const FT &a02, const FT &a03, const FT &a04, const FT &a10, const FT &a11, const FT &a12, const FT &a13, const FT &a14, const FT &a20, const FT &a21, const FT &a22, const FT &a23, const FT &a24, const FT &a30, const FT &a31, const FT &a32, const FT &a33, const FT &a34, const FT &a40, const FT &a41, const FT &a42, const FT &a43, const FT &a44)
 
template<class FT >
FT det6x6_by_formula (const FT &a00, const FT &a01, const FT &a02, const FT &a03, const FT &a04, const FT &a05, const FT &a10, const FT &a11, const FT &a12, const FT &a13, const FT &a14, const FT &a15, const FT &a20, const FT &a21, const FT &a22, const FT &a23, const FT &a24, const FT &a25, const FT &a30, const FT &a31, const FT &a32, const FT &a33, const FT &a34, const FT &a35, const FT &a40, const FT &a41, const FT &a42, const FT &a43, const FT &a44, const FT &a45, const FT &a50, const FT &a51, const FT &a52, const FT &a53, const FT &a54, const FT &a55)
 

Function Documentation

CGAL_BEGIN_NAMESPACE FT det2x2_by_formula ( const FT &  a00,
const FT &  a01,
const FT &  a10,
const FT &  a11 
)
inline

Definition at line 59 of file determinant.h.

Referenced by circumcenter_translateC2(), compare_xC2(), compare_y_at_xC2(), Aff_transformation_repS3< FT >::inverse(), and scaled_distance_to_lineC2().

62 {
63 // First compute the det2x2
64  const FT m01 = a00*a11 - a10*a01;
65  return m01;
66 }

Here is the caller graph for this function:

CGAL_KERNEL_MEDIUM_INLINE FT det3x3_by_formula ( const FT &  a00,
const FT &  a01,
const FT &  a02,
const FT &  a10,
const FT &  a11,
const FT &  a12,
const FT &  a20,
const FT &  a21,
const FT &  a22 
)

Definition at line 71 of file determinant.h.

Referenced by circumcenterC3(), gp_linear_intersection(), Aff_transformation_repS3< FT >::inverse(), scaled_distance_to_planeC3(), and sign_of_determinant3x3().

75 {
76 // First compute the det2x2
77  const FT m01 = a00*a11 - a10*a01;
78  const FT m02 = a00*a21 - a20*a01;
79  const FT m12 = a10*a21 - a20*a11;
80 // Now compute the minors of rank 3
81  const FT m012 = m01*a22 - m02*a12 + m12*a02;
82  return m012;
83 }

Here is the caller graph for this function:

CGAL_KERNEL_LARGE_INLINE FT det4x4_by_formula ( const FT &  a00,
const FT &  a01,
const FT &  a02,
const FT &  a03,
const FT &  a10,
const FT &  a11,
const FT &  a12,
const FT &  a13,
const FT &  a20,
const FT &  a21,
const FT &  a22,
const FT &  a23,
const FT &  a30,
const FT &  a31,
const FT &  a32,
const FT &  a33 
)

Definition at line 88 of file determinant.h.

Referenced by sign_of_determinant4x4().

93 {
94 // First compute the det2x2
95  const FT m01 = a10*a01 - a00*a11;
96  const FT m02 = a20*a01 - a00*a21;
97  const FT m03 = a30*a01 - a00*a31;
98  const FT m12 = a20*a11 - a10*a21;
99  const FT m13 = a30*a11 - a10*a31;
100  const FT m23 = a30*a21 - a20*a31;
101 // Now compute the minors of rank 3
102  const FT m012 = m12*a02 - m02*a12 + m01*a22;
103  const FT m013 = m13*a02 - m03*a12 + m01*a32;
104  const FT m023 = m23*a02 - m03*a22 + m02*a32;
105  const FT m123 = m23*a12 - m13*a22 + m12*a32;
106 // Now compute the minors of rank 4
107  const FT m0123 = m123*a03 - m023*a13 + m013*a23 - m012*a33;
108  return m0123;
109 }

Here is the caller graph for this function:

CGAL_KERNEL_LARGE_INLINE FT det5x5_by_formula ( const FT &  a00,
const FT &  a01,
const FT &  a02,
const FT &  a03,
const FT &  a04,
const FT &  a10,
const FT &  a11,
const FT &  a12,
const FT &  a13,
const FT &  a14,
const FT &  a20,
const FT &  a21,
const FT &  a22,
const FT &  a23,
const FT &  a24,
const FT &  a30,
const FT &  a31,
const FT &  a32,
const FT &  a33,
const FT &  a34,
const FT &  a40,
const FT &  a41,
const FT &  a42,
const FT &  a43,
const FT &  a44 
)

Definition at line 114 of file determinant.h.

Referenced by sign_of_determinant5x5().

120 {
121 // First compute the det2x2
122  const FT m01 = a10*a01 - a00*a11;
123  const FT m02 = a20*a01 - a00*a21;
124  const FT m03 = a30*a01 - a00*a31;
125  const FT m04 = a40*a01 - a00*a41;
126  const FT m12 = a20*a11 - a10*a21;
127  const FT m13 = a30*a11 - a10*a31;
128  const FT m14 = a40*a11 - a10*a41;
129  const FT m23 = a30*a21 - a20*a31;
130  const FT m24 = a40*a21 - a20*a41;
131  const FT m34 = a40*a31 - a30*a41;
132 // Now compute the minors of rank 3
133  const FT m012 = m12*a02 - m02*a12 + m01*a22;
134  const FT m013 = m13*a02 - m03*a12 + m01*a32;
135  const FT m014 = m14*a02 - m04*a12 + m01*a42;
136  const FT m023 = m23*a02 - m03*a22 + m02*a32;
137  const FT m024 = m24*a02 - m04*a22 + m02*a42;
138  const FT m034 = m34*a02 - m04*a32 + m03*a42;
139  const FT m123 = m23*a12 - m13*a22 + m12*a32;
140  const FT m124 = m24*a12 - m14*a22 + m12*a42;
141  const FT m134 = m34*a12 - m14*a32 + m13*a42;
142  const FT m234 = m34*a22 - m24*a32 + m23*a42;
143 // Now compute the minors of rank 4
144  const FT m0123 = m123*a03 - m023*a13 + m013*a23 - m012*a33;
145  const FT m0124 = m124*a03 - m024*a13 + m014*a23 - m012*a43;
146  const FT m0134 = m134*a03 - m034*a13 + m014*a33 - m013*a43;
147  const FT m0234 = m234*a03 - m034*a23 + m024*a33 - m023*a43;
148  const FT m1234 = m234*a13 - m134*a23 + m124*a33 - m123*a43;
149 // Now compute the minors of rank 5
150  const FT m01234 = m1234*a04 - m0234*a14 + m0134*a24 - m0124*a34 + m0123*a44;
151  return m01234;
152 }

Here is the caller graph for this function:

FT det6x6_by_formula ( const FT &  a00,
const FT &  a01,
const FT &  a02,
const FT &  a03,
const FT &  a04,
const FT &  a05,
const FT &  a10,
const FT &  a11,
const FT &  a12,
const FT &  a13,
const FT &  a14,
const FT &  a15,
const FT &  a20,
const FT &  a21,
const FT &  a22,
const FT &  a23,
const FT &  a24,
const FT &  a25,
const FT &  a30,
const FT &  a31,
const FT &  a32,
const FT &  a33,
const FT &  a34,
const FT &  a35,
const FT &  a40,
const FT &  a41,
const FT &  a42,
const FT &  a43,
const FT &  a44,
const FT &  a45,
const FT &  a50,
const FT &  a51,
const FT &  a52,
const FT &  a53,
const FT &  a54,
const FT &  a55 
)

Definition at line 156 of file determinant.h.

Referenced by sign_of_determinant6x6().

163 {
164 // First compute the det2x2
165  const FT m01 = a00*a11 - a10*a01;
166  const FT m02 = a00*a21 - a20*a01;
167  const FT m03 = a00*a31 - a30*a01;
168  const FT m04 = a00*a41 - a40*a01;
169  const FT m05 = a00*a51 - a50*a01;
170  const FT m12 = a10*a21 - a20*a11;
171  const FT m13 = a10*a31 - a30*a11;
172  const FT m14 = a10*a41 - a40*a11;
173  const FT m15 = a10*a51 - a50*a11;
174  const FT m23 = a20*a31 - a30*a21;
175  const FT m24 = a20*a41 - a40*a21;
176  const FT m25 = a20*a51 - a50*a21;
177  const FT m34 = a30*a41 - a40*a31;
178  const FT m35 = a30*a51 - a50*a31;
179  const FT m45 = a40*a51 - a50*a41;
180 // Now compute the minors of rank 3
181  const FT m012 = m01*a22 - m02*a12 + m12*a02;
182  const FT m013 = m01*a32 - m03*a12 + m13*a02;
183  const FT m014 = m01*a42 - m04*a12 + m14*a02;
184  const FT m015 = m01*a52 - m05*a12 + m15*a02;
185  const FT m023 = m02*a32 - m03*a22 + m23*a02;
186  const FT m024 = m02*a42 - m04*a22 + m24*a02;
187  const FT m025 = m02*a52 - m05*a22 + m25*a02;
188  const FT m034 = m03*a42 - m04*a32 + m34*a02;
189  const FT m035 = m03*a52 - m05*a32 + m35*a02;
190  const FT m045 = m04*a52 - m05*a42 + m45*a02;
191  const FT m123 = m12*a32 - m13*a22 + m23*a12;
192  const FT m124 = m12*a42 - m14*a22 + m24*a12;
193  const FT m125 = m12*a52 - m15*a22 + m25*a12;
194  const FT m134 = m13*a42 - m14*a32 + m34*a12;
195  const FT m135 = m13*a52 - m15*a32 + m35*a12;
196  const FT m145 = m14*a52 - m15*a42 + m45*a12;
197  const FT m234 = m23*a42 - m24*a32 + m34*a22;
198  const FT m235 = m23*a52 - m25*a32 + m35*a22;
199  const FT m245 = m24*a52 - m25*a42 + m45*a22;
200  const FT m345 = m34*a52 - m35*a42 + m45*a32;
201 // Now compute the minors of rank 4
202  const FT m0123 = m012*a33 - m013*a23 + m023*a13 - m123*a03;
203  const FT m0124 = m012*a43 - m014*a23 + m024*a13 - m124*a03;
204  const FT m0125 = m012*a53 - m015*a23 + m025*a13 - m125*a03;
205  const FT m0134 = m013*a43 - m014*a33 + m034*a13 - m134*a03;
206  const FT m0135 = m013*a53 - m015*a33 + m035*a13 - m135*a03;
207  const FT m0145 = m014*a53 - m015*a43 + m045*a13 - m145*a03;
208  const FT m0234 = m023*a43 - m024*a33 + m034*a23 - m234*a03;
209  const FT m0235 = m023*a53 - m025*a33 + m035*a23 - m235*a03;
210  const FT m0245 = m024*a53 - m025*a43 + m045*a23 - m245*a03;
211  const FT m0345 = m034*a53 - m035*a43 + m045*a33 - m345*a03;
212  const FT m1234 = m123*a43 - m124*a33 + m134*a23 - m234*a13;
213  const FT m1235 = m123*a53 - m125*a33 + m135*a23 - m235*a13;
214  const FT m1245 = m124*a53 - m125*a43 + m145*a23 - m245*a13;
215  const FT m1345 = m134*a53 - m135*a43 + m145*a33 - m345*a13;
216  const FT m2345 = m234*a53 - m235*a43 + m245*a33 - m345*a23;
217 // Now compute the minors of rank 5
218  const FT m01234 = m0123*a44 - m0124*a34 + m0134*a24 - m0234*a14 + m1234*a04;
219  const FT m01235 = m0123*a54 - m0125*a34 + m0135*a24 - m0235*a14 + m1235*a04;
220  const FT m01245 = m0124*a54 - m0125*a44 + m0145*a24 - m0245*a14 + m1245*a04;
221  const FT m01345 = m0134*a54 - m0135*a44 + m0145*a34 - m0345*a14 + m1345*a04;
222  const FT m02345 = m0234*a54 - m0235*a44 + m0245*a34 - m0345*a24 + m2345*a04;
223  const FT m12345 = m1234*a54 - m1235*a44 + m1245*a34 - m1345*a24 + m2345*a14;
224 // Now compute the minors of rank 6
225  const FT m012345 = m01234*a55 - m01235*a45 + m01245*a35 - m01345*a25 + m02345*a15 - m12345*a05;
226  return m012345;
227 }

Here is the caller graph for this function: