This is the mail archive of the gcc-patches@gcc.gnu.org mailing list for the GCC project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

[PATCH] PR rtl-optimization/37922: Replace read clobbers hard regs.


Replacing a read with shifting operations in dse.c occasionally uses
instructions that set hard regs (like the CC).  This patch makes sure
that the hard regs set by the new code were not already live.


2008-12-17  Kenneth Zadeck <zadeck@naturalbridge.com>

    PR rtl-optimization/37922
    * dce.c (bb_info): Added regs_live field.
    (look_for_hardregs): New function.
    (replace_read): Added regs_live parameter and code to check that
    shift sequence does not clobber live hardregs.
    (check_mem_read_rtx): Added parameter to replace_read.
    (dse_step1): Added regs_live bitmap and initialize it.
    (rest_of_handle_dse): Added DF_NOTES problem and earlier call to
    df_analyze.
    * df-problems.c (df_simulate_one_insn): Renamed to
    df_simulate_one_insn_backwards.
    (df_simulate_artificial_defs_at_top,
    df_simulate_one_insn_forwards,
    df_simulate_artificial_defs_at_end): New functions.
    * df.h (df_simulate_one_insn): Renamed to
    df_simulate_one_insn_backwards.
    (df_simulate_artificial_defs_at_top,
    df_simulate_one_insn_forwards,
    df_simulate_artificial_defs_at_end): New functions.
    * sel-sched.c (propagate_lv_set): Renamed to
    df_simulate_one_insn_backwards.
    * ifcvt.c (dead_or_predicable): Ditto.
    * recog.c (peephole2_optimize): Ditto.
    * rtl-factoring (collect_pattern_seqs, clear_regs_live_in_seq): Ditto.

2008-12-17  Kenneth Zadeck <zadeck@naturalbridge.com>

    PR rtl-optimization/37922
    * g++.dg/torture/pr37922.C: New test.

Bootstrapped and regression tested on x86-32 and x86-64. 
OK for commit?

Index: dse.c
===================================================================
--- dse.c	(revision 142782)
+++ dse.c	(working copy)
@@ -374,6 +374,10 @@ struct bb_info
      operations.  */
   bool apply_wild_read;
 
+  /* The following 4 bitvectors hold information about which positions
+     of which stores are live or dead.  They are indexed by
+     get_bitmap_index.  */
+
   /* The set of store positions that exist in this block before a wild read.  */
   bitmap gen;
   
@@ -401,6 +405,14 @@ struct bb_info
      just initializes the vector from one of the out sets of the
      successors of the block.  */
   bitmap out;
+
+  /* The following bitvector is indexed by the reg number.  It
+     contains the set of regs that are live at the current instruction
+     being processed.  While it contains info for all of the
+     registers, only the pseudos are actually examined.  It is used to
+     assure that shift sequences that are inserted do not accidently
+     clobber live hard regs.  */
+  bitmap regs_live;
 };
 
 typedef struct bb_info *bb_info_t;
@@ -1533,6 +1545,25 @@ find_shift_sequence (int access_size,
 }
 
 
+/* Call back for note_stores to find the hard regs set or clobbered by
+   insn.  Data is a bitmap of the hardregs set so far.  */
+
+static void
+look_for_hardregs (rtx x, const_rtx pat ATTRIBUTE_UNUSED, void *data)
+{
+  bitmap regs_set = (bitmap) data;
+
+  if (REG_P (x)
+      && REGNO (x) < FIRST_PSEUDO_REGISTER)
+    {
+      int regno = REGNO (x);
+      int n = hard_regno_nregs[regno][GET_MODE (x)];
+      while (--n >= 0)
+	bitmap_set_bit (regs_set, regno + n);
+    }
+}
+
+
 /* Take a sequence of:
      A <- r1
      ...
@@ -1566,13 +1597,14 @@ find_shift_sequence (int access_size,
 
 static bool
 replace_read (store_info_t store_info, insn_info_t store_insn, 
-	      read_info_t read_info, insn_info_t read_insn, rtx *loc)
+	      read_info_t read_info, insn_info_t read_insn, rtx *loc, bitmap regs_live)
 {
   enum machine_mode store_mode = GET_MODE (store_info->mem);
   enum machine_mode read_mode = GET_MODE (read_info->mem);
   int shift;
   int access_size; /* In bytes.  */
-  rtx insns, read_reg;
+  rtx insns, this_insn, read_reg;
+  bitmap regs_set;
 
   if (!dbg_cnt (dse))
     return false;
@@ -1624,6 +1656,30 @@ replace_read (store_info_t store_info, i
   insns = get_insns ();
   end_sequence ();
 
+  /* Now we have to scan the set of new instructions to see if the
+     sequence contains and sets of hardregs that happened to be live
+     at this point.  For instance, this can happen if one of the insns
+     sets the CC and the CC happened to be live at that point.  This
+     does occasionally happen, see PR 37922.*/
+
+  regs_set = BITMAP_ALLOC (NULL);
+  for (this_insn = insns; this_insn != NULL_RTX; this_insn = NEXT_INSN (this_insn))
+    note_stores (PATTERN (this_insn), look_for_hardregs, regs_set);
+
+  bitmap_and_into (regs_set, regs_live);
+  if (!bitmap_empty_p (regs_set))
+    {
+      if (dump_file)
+	{
+	  fprintf (dump_file, "abandoning replacement because sequence clobbers live hardregs:");
+	  df_print_regset (dump_file, regs_set);
+	}
+
+      BITMAP_FREE (regs_set);
+      return false;
+    }
+  BITMAP_FREE (regs_set);
+
   if (validate_change (read_insn->insn, loc, read_reg, 0))
     {
       deferred_change_t deferred_change =
@@ -1842,7 +1898,7 @@ check_mem_read_rtx (rtx *loc, void *data
 
 		      if ((store_info->positions_needed & mask) == mask
 			  && replace_read (store_info, i_ptr, 
-					   read_info, insn_info, loc))
+					   read_info, insn_info, loc, bb_info->regs_live))
 			return 0;
 		    }
 		  /* The bases are the same, just see if the offsets
@@ -1911,7 +1967,7 @@ check_mem_read_rtx (rtx *loc, void *data
 
 	      if ((store_info->positions_needed & mask) == mask
 		  && replace_read (store_info, i_ptr, 
-				   read_info, insn_info, loc))
+				   read_info, insn_info, loc, bb_info->regs_live))
 		return 0;
 	    }
 
@@ -2139,7 +2195,8 @@ static void
 dse_step1 (void)
 {
   basic_block bb;
-
+  bitmap regs_live = BITMAP_ALLOC (NULL);
+  
   cselib_init (false);
   all_blocks = BITMAP_ALLOC (NULL);
   bitmap_set_bit (all_blocks, ENTRY_BLOCK);
@@ -2152,6 +2209,10 @@ dse_step1 (void)
 
       memset (bb_info, 0, sizeof (struct bb_info));
       bitmap_set_bit (all_blocks, bb->index);
+      bb_info->regs_live = regs_live;
+
+      bitmap_copy (regs_live, DF_LR_IN (bb));
+      df_simulate_artificial_defs_at_top (bb, regs_live);
 
       bb_table[bb->index] = bb_info;
       cselib_discard_hook = remove_useless_values;
@@ -2172,6 +2233,8 @@ dse_step1 (void)
 	      if (INSN_P (insn))
 		scan_insn (bb_info, insn);
 	      cselib_process_insn (insn);
+	      if (INSN_P (insn))
+		df_simulate_one_insn_forwards (bb, insn, regs_live);
 	    }
 	  
 	  /* This is something of a hack, because the global algorithm
@@ -2238,8 +2301,10 @@ dse_step1 (void)
 
 	  free_alloc_pool (cse_store_info_pool);
 	}
+      bb_info->regs_live = NULL;
     }
 
+  BITMAP_FREE (regs_live);
   cselib_finish ();
   htab_empty (rtx_group_table);
 }
@@ -3239,6 +3304,11 @@ rest_of_handle_dse (void)
 
   df_set_flags (DF_DEFER_INSN_RESCAN);
 
+  /* Need the notes since we must track live hardregs in the forwards
+     direction.  */
+  df_note_add_problem ();
+  df_analyze ();
+
   dse_step0 ();
   dse_step1 ();
   dse_step2_init ();
Index: df-problems.c
===================================================================
--- df-problems.c	(revision 142782)
+++ df-problems.c	(working copy)
@@ -3761,16 +3761,11 @@ df_simulate_fixup_sets (basic_block bb,
 
    df_simulate_artificial_refs_at_end should be called first with a
    bitvector copyied from the DF_LIVE_OUT or DF_LR_OUT.  Then
-   df_simulate_one_insn should be called for each insn in the block,
-   starting with the last on.  Finally,
+   df_simulate_one_insn_backwards should be called for each insn in
+   the block, starting with the last on.  Finally,
    df_simulate_artificial_refs_at_top can be called to get a new value
    of the sets at the top of the block (this is rarely used).
-
-   It would be not be difficult to define a similar set of functions
-   that work in the forwards direction.  In that case the functions
-   would ignore the use sets and look for the REG_DEAD and REG_UNUSED
-   notes.
-----------------------------------------------------------------------------*/
+   ----------------------------------------------------------------------------*/
 
 /* Apply the artificial uses and defs at the end of BB in a backwards
    direction.  */
@@ -3801,7 +3796,7 @@ df_simulate_artificial_refs_at_end (basi
 /* Simulate the backwards effects of INSN on the bitmap LIVE.  */
 
 void 
-df_simulate_one_insn (basic_block bb, rtx insn, bitmap live)
+df_simulate_one_insn_backwards (basic_block bb, rtx insn, bitmap live)
 {
   if (! INSN_P (insn))
     return;	
@@ -3840,3 +3835,89 @@ df_simulate_artificial_refs_at_top (basi
     }
 #endif
 }
+/*----------------------------------------------------------------------------
+   The following three functions are used only for FORWARDS scanning:
+   i.e. they process the defs and the REG_DEAD and REG_UNUSED notes.
+   Thus it is important to add the DF_NOTES problem to the stack of 
+   problems computed before using these functions.
+
+   df_simulate_artificial_defs_at_top should be called first with a
+   bitvector copyied from the DF_LIVE_IN or DF_LR_IN.  Then
+   df_simulate_one_insn_forwards should be called for each insn in
+   the block, starting with the last on.  Finally,
+   df_simulate_artificial_def_at_end can be called to get a new value
+   of the sets at the top of the block (this is rarely used).
+   ----------------------------------------------------------------------------*/
+
+/* Apply the artificial uses and defs at the top of BB in a backwards
+   direction.  */
+
+void 
+df_simulate_artificial_defs_at_top (basic_block bb, bitmap live)
+{
+  df_ref *def_rec;
+  int bb_index = bb->index;
+  
+  for (def_rec = df_get_artificial_defs (bb_index); *def_rec; def_rec++)
+    {
+      df_ref def = *def_rec;
+      if (DF_REF_FLAGS (def) & DF_REF_AT_TOP)
+	bitmap_clear_bit (live, DF_REF_REGNO (def));
+    }
+}
+
+/* Simulate the backwards effects of INSN on the bitmap LIVE.  */
+
+void 
+df_simulate_one_insn_forwards (basic_block bb, rtx insn, bitmap live)
+{
+  rtx link;
+  if (! INSN_P (insn))
+    return;	
+  
+  df_simulate_defs (insn, live);
+
+  /* Clear all of the registers that go dead.  */
+  for (link = REG_NOTES (insn); link; link = XEXP (link, 1))
+    {
+      switch (REG_NOTE_KIND (link))
+	case REG_DEAD:
+	case REG_UNUSED:
+	{
+	  rtx reg = XEXP (link, 0);
+	  int regno = REGNO (reg);
+	  if (regno < FIRST_PSEUDO_REGISTER)
+	    {
+	      int n = hard_regno_nregs[regno][GET_MODE (reg)];
+	      while (--n >= 0)
+		bitmap_clear_bit (live, regno + n);
+	    }
+	  else 
+	    bitmap_clear_bit (live, regno);
+	  break;
+	default:
+	  break;
+	}
+    }
+  df_simulate_fixup_sets (bb, live);
+}
+
+
+/* Apply the artificial uses and defs at the end of BB in a backwards
+   direction.  */
+
+void 
+df_simulate_artificial_defs_at_end (basic_block bb, bitmap live)
+{
+  df_ref *def_rec;
+  int bb_index = bb->index;
+  
+  for (def_rec = df_get_artificial_defs (bb_index); *def_rec; def_rec++)
+    {
+      df_ref def = *def_rec;
+      if ((DF_REF_FLAGS (def) & DF_REF_AT_TOP) == 0)
+	bitmap_clear_bit (live, DF_REF_REGNO (def));
+    }
+}
+
+
Index: df.h
===================================================================
--- df.h	(revision 142782)
+++ df.h	(working copy)
@@ -960,8 +960,11 @@ extern void df_simulate_find_defs (rtx,
 extern void df_simulate_defs (rtx, bitmap);
 extern void df_simulate_uses (rtx, bitmap);
 extern void df_simulate_artificial_refs_at_end (basic_block, bitmap);
-extern void df_simulate_one_insn (basic_block, rtx, bitmap);
+extern void df_simulate_one_insn_backwards (basic_block, rtx, bitmap);
 extern void df_simulate_artificial_refs_at_top (basic_block, bitmap);
+extern void df_simulate_artificial_defs_at_end (basic_block, bitmap);
+extern void df_simulate_one_insn_forwards (basic_block, rtx, bitmap);
+extern void df_simulate_artificial_defs_at_top (basic_block, bitmap);
 
 /* Functions defined in df-scan.c.  */
 
Index: sel-sched.c
===================================================================
--- sel-sched.c	(revision 142782)
+++ sel-sched.c	(working copy)
@@ -2947,7 +2947,7 @@ propagate_lv_set (regset lv, insn_t insn
   if (INSN_NOP_P (insn))
     return;
 
-  df_simulate_one_insn (BLOCK_FOR_INSN (insn), insn, lv);
+  df_simulate_one_insn_backwards (BLOCK_FOR_INSN (insn), insn, lv);
 }
 
 /* Return livness set at the end of BB.  */
Index: ifcvt.c
===================================================================
--- ifcvt.c	(revision 142782)
+++ ifcvt.c	(working copy)
@@ -3956,7 +3956,7 @@ dead_or_predicable (basic_block test_bb,
 	  if (INSN_P (insn))
 	    {
 	      df_simulate_find_defs (insn, test_set);
-	      df_simulate_one_insn (test_bb, insn, test_live);
+	      df_simulate_one_insn_backwards (test_bb, insn, test_live);
 	    }
 	  prev = PREV_INSN (insn);
 	  if (insn == earliest)
Index: recog.c
===================================================================
--- recog.c	(revision 142782)
+++ recog.c	(working copy)
@@ -3059,7 +3059,7 @@ peephole2_optimize (void)
 		  && peep2_insn_data[peep2_current].insn == NULL_RTX)
 		peep2_current_count++;
 	      peep2_insn_data[peep2_current].insn = insn;
-	      df_simulate_one_insn (bb, insn, live);
+	      df_simulate_one_insn_backwards (bb, insn, live);
 	      COPY_REG_SET (peep2_insn_data[peep2_current].live_before, live);
 
 	      if (RTX_FRAME_RELATED_P (insn))
@@ -3218,7 +3218,7 @@ peephole2_optimize (void)
 			    peep2_current_count++;
 			  peep2_insn_data[i].insn = x;
 			  df_insn_rescan (x);
-			  df_simulate_one_insn (bb, x, live);
+			  df_simulate_one_insn_backwards (bb, x, live);
 			  bitmap_copy (peep2_insn_data[i].live_before, live);
 			}
 		      x = PREV_INSN (x);
Index: rtl-factoring.c
===================================================================
--- rtl-factoring.c	(revision 142782)
+++ rtl-factoring.c	(working copy)
@@ -486,7 +486,7 @@ collect_pattern_seqs (void)
 	  }
 	if (insn == BB_HEAD (bb))
 	  break;
-	df_simulate_one_insn (bb, insn, &live);
+	df_simulate_one_insn_backwards (bb, insn, &live);
 	insn = prev;
       }
 
@@ -576,7 +576,7 @@ clear_regs_live_in_seq (HARD_REG_SET * r
 
   /* Propagate until INSN if found.  */
   for (x = BB_END (bb); x != insn; x = PREV_INSN (x))
-    df_simulate_one_insn (bb, x, &live);
+    df_simulate_one_insn_backwards (bb, x, &live);
 
   /* Clear registers live after INSN.  */
   renumbered_reg_set_to_hard_reg_set (&hlive, &live);
@@ -586,7 +586,7 @@ clear_regs_live_in_seq (HARD_REG_SET * r
   for (i = 0; i < length;)
     {
       rtx prev = PREV_INSN (x);
-      df_simulate_one_insn (bb, x, &live);
+      df_simulate_one_insn_backwards (bb, x, &live);
 
       if (INSN_P (x))
         {
Index: testsuite/g++.dg/torture/pr37922.C
===================================================================
--- testsuite/g++.dg/torture/pr37922.C	(revision 0)
+++ testsuite/g++.dg/torture/pr37922.C	(revision 0)
@@ -0,0 +1,502 @@
+// { dg-do run }
+// { dg-options "-fpic" { target fpic } }
+
+typedef __SIZE_TYPE__ size_t;
+
+template <typename NumType>
+inline
+NumType
+absolute(NumType const& x)
+{
+  if (x < NumType(0)) return -x;
+  return x;
+}
+
+class trivial_accessor
+{
+  public:
+    typedef size_t index_type;
+    struct index_value_type {};
+
+    trivial_accessor() : size_(0) {}
+
+    trivial_accessor(size_t const& n) : size_(n) {}
+
+    size_t size_1d() const { return size_; }
+
+  protected:
+    size_t size_;
+};
+
+namespace N0
+{
+  template <typename ElementType,
+            typename AccessorType = trivial_accessor>
+  class const_ref
+  {
+    public:
+      typedef ElementType value_type;
+      typedef size_t size_type;
+
+      typedef AccessorType accessor_type;
+      typedef typename accessor_type::index_type index_type;
+      typedef typename accessor_type::index_value_type index_value_type;
+
+      const_ref() {}
+
+      const_ref(const ElementType* begin, accessor_type const& accessor)
+      : begin_(begin), accessor_(accessor)
+      {
+        init();
+      }
+
+      const_ref(const ElementType* begin, index_value_type const& n0)
+      : begin_(begin), accessor_(n0)
+      {
+        init();
+      }
+
+      const_ref(const ElementType* begin, index_value_type const& n0,
+                                          index_value_type const& n1)
+      : begin_(begin), accessor_(n0, n1)
+      {
+        init();
+      }
+
+      const_ref(const ElementType* begin, index_value_type const& n0,
+                                          index_value_type const& n1,
+                                          index_value_type const& n2)
+      : begin_(begin), accessor_(n0, n1, n2)
+      {
+        init();
+      }
+
+      accessor_type const& accessor() const { return accessor_; }
+      size_type size() const { return size_; }
+
+      const ElementType* begin() const { return begin_; }
+      const ElementType* end() const { return end_; }
+
+      ElementType const&
+      operator[](size_type i) const { return begin_[i]; }
+
+      const_ref<ElementType>
+      as_1d() const
+      {
+        return const_ref<ElementType>(begin_, size_);
+      }
+
+    protected:
+      void
+      init()
+      {
+        size_ = accessor_.size_1d();
+        end_ = begin_ + size_;
+      }
+
+      const ElementType* begin_;
+      accessor_type accessor_;
+      size_type size_;
+      const ElementType* end_;
+  };
+}
+
+template <typename ElementType,
+          typename AccessorType = trivial_accessor>
+class ref : public N0::const_ref<ElementType, AccessorType>
+{
+  public:
+    typedef ElementType value_type;
+    typedef size_t size_type;
+
+    typedef N0::const_ref<ElementType, AccessorType> base_class;
+    typedef AccessorType accessor_type;
+    typedef typename accessor_type::index_type index_type;
+
+    ref() {}
+
+    ElementType*
+    begin() const { return const_cast<ElementType*>(this->begin_); }
+
+    ElementType*
+    end() const { return const_cast<ElementType*>(this->end_); }
+
+    ElementType&
+    operator[](size_type i) const { return begin()[i]; }
+};
+
+namespace N1 {
+  template <typename ElementType, size_t N>
+  class tiny_plain
+  {
+    public:
+      typedef ElementType value_type;
+      typedef size_t size_type;
+
+      static const size_t fixed_size=N;
+
+      ElementType elems[N];
+
+      tiny_plain() {}
+
+      static size_type size() { return N; }
+
+      ElementType* begin() { return elems; }
+      const ElementType* begin() const { return elems; }
+      ElementType* end() { return elems+N; }
+      const ElementType* end() const { return elems+N; }
+      ElementType& operator[](size_type i) { return elems[i]; }
+      ElementType const& operator[](size_type i) const { return elems[i]; }
+  };
+
+  template <typename ElementType, size_t N>
+  class tiny : public tiny_plain<ElementType, N>
+  {
+    public:
+      typedef ElementType value_type;
+      typedef size_t size_type;
+
+      typedef tiny_plain<ElementType, N> base_class;
+
+      tiny() {}
+  };
+}
+
+template <typename NumType>
+class mat3 : public N1::tiny_plain<NumType, 9>
+{
+  public:
+    typedef typename N1::tiny_plain<NumType, 9> base_type;
+
+    mat3() {}
+    mat3(NumType const& e00, NumType const& e01, NumType const& e02,
+         NumType const& e10, NumType const& e11, NumType const& e12,
+         NumType const& e20, NumType const& e21, NumType const& e22)
+      : base_type(e00, e01, e02, e10, e11, e12, e20, e21, e22)
+    {}
+    mat3(base_type const& a)
+      : base_type(a)
+    {}
+
+    NumType const&
+    operator()(size_t r, size_t c) const
+    {
+      return this->elems[r * 3 + c];
+    }
+    NumType&
+    operator()(size_t r, size_t c)
+    {
+      return this->elems[r * 3 + c];
+    }
+
+    NumType
+    trace() const
+    {
+      mat3 const& m = *this;
+      return m[0] + m[4] + m[8];
+    }
+
+    NumType
+    determinant() const
+    {
+      mat3 const& m = *this;
+      return   m[0] * (m[4] * m[8] - m[5] * m[7])
+             - m[1] * (m[3] * m[8] - m[5] * m[6])
+             + m[2] * (m[3] * m[7] - m[4] * m[6]);
+    }
+};
+
+template <typename NumType>
+inline
+mat3<NumType>
+operator-(mat3<NumType> const& v)
+{
+  mat3<NumType> result;
+  for(size_t i=0;i<9;i++) {
+    result[i] = -v[i];
+  }
+  return result;
+}
+
+class mat_grid : public N1::tiny<size_t, 2>
+{
+  public:
+    typedef N1::tiny<size_t, 2> index_type;
+    typedef index_type::value_type index_value_type;
+
+    mat_grid() { this->elems[0]=0; this->elems[1]=0; }
+
+    mat_grid(index_type const& n) : index_type(n) {}
+
+    mat_grid(index_value_type const& n0, index_value_type const& n1)
+    { this->elems[0]=n0; this->elems[1]=n1; }
+
+    size_t size_1d() const { return elems[0] * elems[1]; }
+
+    size_t
+    operator()(index_value_type const& r, index_value_type const& c) const
+    {
+      return r * elems[1] + c;
+    }
+};
+
+template <typename NumType, typename AccessorType = mat_grid>
+class mat_const_ref : public N0::const_ref<NumType, AccessorType>
+{
+  public:
+    typedef AccessorType accessor_type;
+    typedef typename N0::const_ref<NumType, AccessorType> base_type;
+    typedef typename accessor_type::index_value_type index_value_type;
+
+    mat_const_ref() {}
+
+    mat_const_ref(const NumType* begin, accessor_type const& grid)
+    : base_type(begin, grid)
+    {}
+
+    mat_const_ref(const NumType* begin, index_value_type const& n_rows,
+                                        index_value_type const& n_columns)
+    : base_type(begin, accessor_type(n_rows, n_columns))
+    {}
+
+    accessor_type
+    grid() const { return this->accessor(); }
+
+    index_value_type const&
+    n_rows() const { return this->accessor()[0]; }
+
+    index_value_type const&
+    n_columns() const { return this->accessor()[1]; }
+
+    NumType const&
+    operator()(index_value_type const& r, index_value_type const& c) const
+    {
+      return this->begin()[this->accessor()(r, c)];
+    }
+};
+
+template <typename NumType, typename AccessorType = mat_grid>
+class mat_ref : public mat_const_ref<NumType, AccessorType>
+{
+  public:
+    typedef AccessorType accessor_type;
+    typedef mat_const_ref<NumType, AccessorType> base_type;
+    typedef typename accessor_type::index_value_type index_value_type;
+
+    mat_ref() {}
+
+    mat_ref(NumType* begin, accessor_type const& grid)
+    : base_type(begin, grid)
+    {}
+
+    mat_ref(NumType* begin, index_value_type n_rows,
+                            index_value_type n_columns)
+    : base_type(begin, accessor_type(n_rows, n_columns))
+    {}
+
+    NumType*
+    begin() const { return const_cast<NumType*>(this->begin_); }
+
+    NumType*
+    end() const { return const_cast<NumType*>(this->end_); }
+
+    NumType&
+    operator[](index_value_type const& i) const { return begin()[i]; }
+
+    NumType&
+    operator()(index_value_type const& r, index_value_type const& c) const
+    {
+      return this->begin()[this->accessor()(r, c)];
+    }
+};
+
+  template <typename AnyType>
+  inline void
+  swap(AnyType* a, AnyType* b, size_t n)
+  {
+    for(size_t i=0;i<n;i++) {
+      AnyType t = a[i]; a[i] = b[i]; b[i] = t;
+    }
+  }
+
+template <typename IntType>
+size_t
+form_t(mat_ref<IntType>& m,
+       mat_ref<IntType> const& t)
+{
+  typedef size_t size_t;
+  size_t mr = m.n_rows();
+  size_t mc = m.n_columns();
+  size_t tc = t.n_columns();
+  if (tc) {
+  }
+  size_t i, j;
+  for (i = j = 0; i < mr && j < mc;) {
+    size_t k = i; while (k < mr && m(k,j) == 0) k++;
+    if (k == mr)
+      j++;
+    else {
+      if (i != k) {
+                swap(&m(i,0), &m(k,0), mc);
+        if (tc) swap(&t(i,0), &t(k,0), tc);
+      }
+      for (k++; k < mr; k++) {
+        IntType a = absolute(m(k, j));
+        if (a != 0 && a < absolute(m(i,j))) {
+                  swap(&m(i,0), &m(k,0), mc);
+          if (tc) swap(&t(i,0), &t(k,0), tc);
+        }
+      }
+      if (m(i,j) < 0) {
+                for(size_t ic=0;ic<mc;ic++) m(i,ic) *= -1;
+        if (tc) for(size_t ic=0;ic<tc;ic++) t(i,ic) *= -1;
+      }
+      bool cleared = true;
+      for (k = i+1; k < mr; k++) {
+        IntType a = m(k,j) / m(i,j);
+        if (a != 0) {
+                  for(size_t ic=0;ic<mc;ic++) m(k,ic) -= a * m(i,ic);
+          if (tc) for(size_t ic=0;ic<tc;ic++) t(k,ic) -= a * t(i,ic);
+        }
+        if (m(k,j) != 0) cleared = false;
+      }
+      if (cleared) { i++; j++; }
+    }
+  }
+  m = mat_ref<IntType>(m.begin(), i, mc);
+  return i;
+}
+
+template <typename IntType>
+size_t
+form(mat_ref<IntType>& m)
+{
+  mat_ref<IntType> t(0,0,0);
+  return form_t(m, t);
+}
+
+typedef mat3<int> sg_mat3;
+
+class rot_mx
+{
+  public:
+    explicit
+    rot_mx(sg_mat3 const& m, int denominator=1)
+    : num_(m), den_(denominator)
+    {}
+
+    sg_mat3 const&
+    num() const { return num_; }
+    sg_mat3&
+    num()       { return num_; }
+
+    int const&
+    operator[](size_t i) const { return num_[i]; }
+    int&
+    operator[](size_t i)       { return num_[i]; }
+
+    int
+    const& operator()(int r, int c) const { return num_(r, c); }
+    int&
+    operator()(int r, int c)       { return num_(r, c); }
+
+    int const&
+    den() const { return den_; }
+    int&
+    den()       { return den_; }
+
+    rot_mx
+    minus_unit_mx() const
+    {
+      rot_mx result(*this);
+      for (size_t i=0;i<9;i+=4) result[i] -= den_;
+      return result;
+    }
+
+    rot_mx
+    operator-() const { return rot_mx(-num_, den_); }
+
+    int
+    type() const;
+
+    int
+    order(int type=0) const;
+
+  private:
+    sg_mat3 num_;
+    int den_;
+};
+
+class rot_mx_info
+{
+  public:
+    rot_mx_info(rot_mx const& r);
+
+    int type() const { return type_; }
+
+  private:
+    int type_;
+};
+
+int rot_mx::type() const
+{
+  int det = num_.determinant();
+  if (det == -1 || det == 1) {
+    switch (num_.trace()) {
+      case -3:                return -1;
+      case -2:                return -6;
+      case -1: if (det == -1) return -4;
+               else           return  2;
+      case  0: if (det == -1) return -3;
+               else           return  3;
+      case  1: if (det == -1) return -2;
+               else           return  4;
+      case  2:                return  6;
+      case  3:                return  1;
+    }
+  }
+  return 0;
+}
+
+int rot_mx::order(int type) const
+{
+  if (type == 0) type = rot_mx::type();
+  if (type > 0) return  type;
+  if (type % 2) return -type * 2;
+                return -type;
+}
+
+rot_mx_info::rot_mx_info(rot_mx const& r)
+: type_(r.type())
+{
+  if (type_ == 0) {
+    return;
+  }
+  rot_mx proper_r = r;
+  int proper_order = type_;
+  // THE PROBLEM IS AROUND HERE
+  if (proper_order < 0) {
+    proper_order *= -1;
+    proper_r = -proper_r; // THIS FAILS ...
+  }
+  if (proper_order > 1) {
+    rot_mx rmi = proper_r.minus_unit_mx(); // ... THEREFORE WRONG HERE
+    mat_ref<int> re_mx(rmi.num().begin(), 3, 3);
+    if (form(re_mx) != 2) {
+      type_ = 0;
+    }
+  }
+}
+
+int main()
+{
+  N1::tiny<int, 9> e;
+  e[0] = 1; e[1] =  0; e[2] = 0;
+  e[3] = 0; e[4] = -1; e[5] = 0;
+  e[6] = 0; e[7] =  0; e[8] = 1;
+  rot_mx r(e);
+  rot_mx_info ri(r);
+  if (ri.type() != -2)
+    __builtin_abort ();
+  return 0;
+}

Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]