Return to hamdist.c CVS log

Up to [local] / OpenXM_contrib / gmp / mpz

Annotation of OpenXM_contrib/gmp/mpz/hamdist.c, Revision 1.1.1.3

1.1.1.3 ! ohara       1: /* mpz_hamdist -- calculate hamming distance.
1.1       maekawa     2:
1.1.1.3 ! ohara       3: Copyright 1994, 1996, 2001, 2002 Free Software Foundation, Inc.
1.1       maekawa     4:
                      5: This file is part of the GNU MP Library.
                      6:
                      7: The GNU MP Library is free software; you can redistribute it and/or modify
1.1.1.2   maekawa     8: it under the terms of the GNU Lesser General Public License as published by
                      9: the Free Software Foundation; either version 2.1 of the License, or (at your
1.1       maekawa    10: option) any later version.
                     11:
                     12: The GNU MP Library is distributed in the hope that it will be useful, but
                     13: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
1.1.1.2   maekawa    14: or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
1.1       maekawa    15: License for more details.
                     16:
1.1.1.2   maekawa    17: You should have received a copy of the GNU Lesser General Public License
1.1       maekawa    18: along with the GNU MP Library; see the file COPYING.LIB.  If not, write to
                     19: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
                     20: MA 02111-1307, USA. */
                     21:
                     22: #include "gmp.h"
                     23: #include "gmp-impl.h"
                     24:
1.1.1.3 ! ohara      25:
        !            26: unsigned long
1.1       maekawa    27: mpz_hamdist (mpz_srcptr u, mpz_srcptr v)
                     28: {
1.1.1.3 ! ohara      29:   mp_srcptr      up, vp;
        !            30:   mp_size_t      usize, vsize;
        !            31:   unsigned long  count;
1.1       maekawa    32:
1.1.1.3 ! ohara      33:   usize = SIZ(u);
        !            34:   vsize = SIZ(v);
1.1       maekawa    35:
1.1.1.3 ! ohara      36:   up = PTR(u);
        !            37:   vp = PTR(v);
1.1       maekawa    38:
1.1.1.3 ! ohara      39:   if (usize >= 0)
1.1       maekawa    40:     {
1.1.1.3 ! ohara      41:       if (vsize < 0)
        !            42:         return ~ (unsigned long) 0;
        !            43:
        !            44:       /* positive/positive */
        !            45:
        !            46:       if (usize < vsize)
        !            47:         MPN_SRCPTR_SWAP (up,usize, vp,vsize);
        !            48:
        !            49:       count = 0;
        !            50:       if (vsize != 0)
        !            51:         count = mpn_hamdist (up, vp, vsize);
        !            52:
        !            53:       usize -= vsize;
        !            54:       if (usize != 0)
        !            55:         count += mpn_popcount (up + vsize, usize);
        !            56:
        !            57:       return count;
1.1       maekawa    58:     }
                     59:   else
                     60:     {
1.1.1.3 ! ohara      61:       mp_limb_t  ulimb, vlimb;
        !            62:       mp_size_t  old_vsize, step;
        !            63:
        !            64:       if (vsize >= 0)
        !            65:         return ~ (unsigned long) 0;
1.1       maekawa    66:
1.1.1.3 ! ohara      67:       /* negative/negative */
        !            68:
        !            69:       usize = -usize;
        !            70:       vsize = -vsize;
        !            71:
        !            72:       /* skip common low zeros */
        !            73:       for (;;)
        !            74:         {
        !            75:           ASSERT (usize > 0);
        !            76:           ASSERT (vsize > 0);
        !            77:
        !            78:           usize--;
        !            79:           vsize--;
        !            80:
        !            81:           ulimb = *up++;
        !            82:           vlimb = *vp++;
        !            83:
        !            84:           if (ulimb != 0)
        !            85:             break;
        !            86:
        !            87:           if (vlimb != 0)
        !            88:             {
        !            89:               MPN_SRCPTR_SWAP (up,usize, vp,vsize);
        !            90:               ulimb = vlimb;
        !            91:               vlimb = 0;
        !            92:               break;
        !            93:             }
        !            94:         }
        !            95:
        !            96:       /* twos complement first non-zero limbs (ulimb is non-zero, but vlimb
        !            97:          might be zero) */
        !            98:       ulimb = -ulimb;
        !            99:       vlimb = -vlimb;
        !           100:       popc_limb (count, (ulimb ^ vlimb) & GMP_NUMB_MASK);
        !           101:
        !           102:       if (vlimb == 0)
        !           103:         {
        !           104:           unsigned long  twoscount;
        !           105:
        !           106:           /* first non-zero of v */
        !           107:           old_vsize = vsize;
        !           108:           do
        !           109:             {
        !           110:               ASSERT (vsize > 0);
        !           111:               vsize--;
        !           112:               vlimb = *vp++;
        !           113:             }
        !           114:           while (vlimb == 0);
        !           115:
        !           116:           /* part of u corresponding to skipped v zeros */
        !           117:           step = old_vsize - vsize - 1;
        !           118:           count += step * GMP_NUMB_BITS;
        !           119:           step = MIN (step, usize);
        !           120:           if (step != 0)
        !           121:             {
        !           122:               count -= mpn_popcount (up, step);
        !           123:               usize -= step;
        !           124:               up += step;
        !           125:             }
        !           126:
        !           127:           /* First non-zero vlimb as twos complement, xor with ones
        !           128:              complement ulimb.  Note -v^(~0^u) == (v-1)^u. */
        !           129:           vlimb--;
        !           130:           if (usize != 0)
        !           131:             {
        !           132:               usize--;
        !           133:               vlimb ^= *up++;
        !           134:             }
        !           135:           popc_limb (twoscount, vlimb);
        !           136:           count += twoscount;
        !           137:         }
        !           138:
        !           139:       /* Overlapping part of u and v, if any.  Ones complement both, so just
        !           140:          plain hamdist. */
        !           141:       step = MIN (usize, vsize);
        !           142:       if (step != 0)
        !           143:         {
        !           144:           count += mpn_hamdist (up, vp, step);
        !           145:           usize -= step;
        !           146:           vsize -= step;
        !           147:           up += step;
        !           148:           vp += step;
        !           149:         }
        !           150:
        !           151:       /* Remaining high part of u or v, if any, ones complement but xor
        !           152:          against all ones in the other, so plain popcount. */
        !           153:       if (usize != 0)
        !           154:         {
        !           155:         remaining:
        !           156:           count += mpn_popcount (up, usize);
        !           157:         }
        !           158:       else if (vsize != 0)
        !           159:         {
        !           160:           up = vp;
        !           161:           usize = vsize;
        !           162:           goto remaining;
        !           163:         }
        !           164:       return count;
        !           165:     }
1.1       maekawa   166: }