]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
constructor for mixed elements
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 15 Feb 1999 16:23:07 +0000 (16:23 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 15 Feb 1999 16:23:07 +0000 (16:23 +0000)
git-svn-id: https://svn.dealii.org/trunk@802 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_system.h

index df5fa8cbced6c7cc105aadff9dc3a798133d0dbd..5a5beaabf80bbbd48a681e164bf9f32d657f8e6d 100644 (file)
@@ -63,7 +63,7 @@
  * function) to the respective component of the finite element, without interaction
  * of the different components.
  *
- * @author Wolfgang Bangerth, 1999
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999
  */
 template <int dim>
 class FESystem : public FiniteElement<dim>
@@ -95,9 +95,20 @@ class FESystem : public FiniteElement<dim>
                                      * class needs to be of the same dimension
                                      * as is this object.
                                      */
-    template <typename FE>
+    template <class FE>
     FESystem (const FE &fe, const unsigned int n_elements);
 
+                                    /** 
+                                     * Constructor for mixed
+                                     * discretizations with two
+                                     * components.
+                                     *
+                                     * See the other constructor.
+                                     */
+    template <class FE1, class FE2>
+    FESystem (const FE1 &fe1, const unsigned int n1,
+             const FE2 &fe2, const unsigned int n2);
+
                                     /**
                                      * Destructor.
                                      */
@@ -377,6 +388,15 @@ class FESystem : public FiniteElement<dim>
     multiply_dof_numbers (const FiniteElementData<dim> &fe_data,
                          const unsigned int            N);
     
+                                    /**
+                                     * Same as above for mixed elements.
+                                     */
+    static FiniteElementData<dim>
+    multiply_dof_numbers (const FiniteElementData<dim> &fe1,
+                         const unsigned int            N1,
+                         const FiniteElementData<dim> &fe2,
+                         const unsigned int            N2);
+    
                                     /**
                                      * This function is simply singled out of
                                      * the constructor. It sets up the
@@ -429,8 +449,23 @@ FESystem<dim>::FESystem (const FE &fe, const unsigned int n_elements) :
                FiniteElement (multiply_dof_numbers(fe, n_elements)),
                base_elements(1)
 {
-  base_elements[0] = ElementPair(&fe, n_elements);
+  base_elements[0] = ElementPair(new FE, n_elements);
+  base_elements[0].first -> subscribe ();
+  initialize ();
+};
+
+template <int dim>
+template <class FE1, class FE2>
+FESystem<dim>::FESystem (const FE1 &fe1, const unsigned int n1,
+                        const FE2 &fe2, const unsigned int n2)
+               :
+               FiniteElement (multiply_dof_numbers(fe1, n1, fe2, n2)),
+               base_elements(2)
+{
+  base_elements[0] = ElementPair(new FE1, n1);
   base_elements[0].first -> subscribe ();
+  base_elements[1] = ElementPair(new FE2, n2);
+  base_elements[1].first -> subscribe ();
   initialize ();
 };
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.