ENH: Handle NaT in most operations, some progress towards proper casting rules

8 files changed, 407 insertions, 62 deletions

diff --git a/numpy/core/code_generators/generate_umath.py b/numpy/core/code_generators/generate_umath.py
index db29b467e..945816ad6 100644
--- a/numpy/core/code_generators/generate_umath.py
+++ b/numpy/core/code_generators/generate_umath.py
@@ -476,7 +476,6 @@ defdict = {
           None,
           TD(flts, f="logaddexp2", astype={'e':'f'})
           ),
-# FIXME: decide if the times should have the bitwise operations.
 'bitwise_and' :
     Ufunc(2, 1, One,
           docstrings.get('numpy.core.umath.bitwise_and'),
diff --git a/numpy/core/src/multiarray/_datetime.h b/numpy/core/src/multiarray/_datetime.h
index 0ad7b3404..dbff004b9 100644
--- a/numpy/core/src/multiarray/_datetime.h
+++ b/numpy/core/src/multiarray/_datetime.h
@@ -93,6 +93,17 @@ parse_datetime_unit_from_string(char *str, Py_ssize_t len, char *metastr);
 NPY_NO_EXPORT int
 convert_datetime_divisor_to_multiple(PyArray_DatetimeMetaData *meta,
                                     int den, char *metastr);
+
+/*
+ * Determines whether the 'divisor' metadata divides evenly into
+ * the 'dividend' metadata.
+ */
+NPY_NO_EXPORT npy_bool
+datetime_metadata_divides(
+                        PyArray_Descr *dividend,
+                        PyArray_Descr *divisor,
+                        int strict_with_nonlinear_units);
+
 /*
  * Computes the GCD of the two date-time metadata values. Raises
  * an exception if there is no reasonable GCD, such as with
diff --git a/numpy/core/src/multiarray/convert_datatype.c b/numpy/core/src/multiarray/convert_datatype.c
index e1fe71c6b..a0d1c76d3 100644
--- a/numpy/core/src/multiarray/convert_datatype.c
+++ b/numpy/core/src/multiarray/convert_datatype.c
@@ -235,20 +235,30 @@ PyArray_CanCastTo(PyArray_Descr *from, PyArray_Descr *to)
     ret = (npy_bool) PyArray_CanCastSafely(fromtype, totype);
     if (ret) {
         /* Check String and Unicode more closely */
-        if (fromtype == PyArray_STRING) {
-            if (totype == PyArray_STRING) {
+        if (fromtype == NPY_STRING) {
+            if (totype == NPY_STRING) {
                 ret = (from->elsize <= to->elsize);
             }
-            else if (totype == PyArray_UNICODE) {
+            else if (totype == NPY_UNICODE) {
                 ret = (from->elsize << 2 <= to->elsize);
             }
         }
-        else if (fromtype == PyArray_UNICODE) {
-            if (totype == PyArray_UNICODE) {
+        else if (fromtype == NPY_UNICODE) {
+            if (totype == NPY_UNICODE) {
                 ret = (from->elsize <= to->elsize);
             }
         }
         /*
+         * For datetime/timedelta, only treat casts moving towards
+         * more precision as safe.
+         */
+        else if (fromtype == NPY_DATETIME && totype == NPY_DATETIME) {
+            return datetime_metadata_divides(from, to, 0);
+        }
+        else if (fromtype == NPY_TIMEDELTA && totype == NPY_TIMEDELTA) {
+            return datetime_metadata_divides(from, to, 1);
+        }
+        /*
          * TODO: If totype is STRING or unicode
          * see if the length is long enough to hold the
          * stringified value of the object.
@@ -290,9 +300,12 @@ dtype_kind_to_ordering(char kind)
         /* Object kind */
         case 'O':
             return 9;
-        /* Anything else - ideally shouldn't happen... */
+        /*
+         * Anything else, like datetime, is special cased to
+         * not fit in this hierarchy
+         */
         default:
-            return 10;
+            return -1;
     }
 }
 
@@ -362,9 +375,7 @@ PyArray_CanCastTypeTo(PyArray_Descr *from, PyArray_Descr *to,
 
         switch (casting) {
             case NPY_NO_CASTING:
-                return (from->elsize == to->elsize) &&
-                        PyArray_ISNBO(from->byteorder) ==
-                                    PyArray_ISNBO(to->byteorder);
+                return PyArray_EquivTypes(from, to);
             case NPY_EQUIV_CASTING:
                 return (from->elsize == to->elsize);
             case NPY_SAFE_CASTING:
@@ -382,9 +393,15 @@ PyArray_CanCastTypeTo(PyArray_Descr *from, PyArray_Descr *to,
             /*
              * Also allow casting from lower to higher kinds, according
              * to the ordering provided by dtype_kind_to_ordering.
+             * Some kinds, like datetime, don't fit in the hierarchy,
+             * and are special cased as -1.
              */
-            return dtype_kind_to_ordering(from->kind) <=
-                            dtype_kind_to_ordering(to->kind);
+            int from_order, to_order;
+
+            from_order = dtype_kind_to_ordering(from->kind);
+            to_order = dtype_kind_to_ordering(to->kind);
+
+            return from_order != -1 && from_order <= to_order;
         }
         else {
             return 0;
@@ -651,10 +668,9 @@ PyArray_PromoteTypes(PyArray_Descr *type1, PyArray_Descr *type2)
     }
 
     switch (type_num1) {
-        /* BOOL can convert to anything except datetime/timedelta/void */
+        /* BOOL can convert to anything except datetime/void */
         case NPY_BOOL:
-            if (type_num2 != NPY_DATETIME && type_num2 != NPY_TIMEDELTA &&
-                                    type_num2 != NPY_VOID) {
+            if (type_num2 != NPY_DATETIME && type_num2 != NPY_VOID) {
                 Py_INCREF(type2);
                 return type2;
             }
diff --git a/numpy/core/src/multiarray/datetime.c b/numpy/core/src/multiarray/datetime.c
index 729cc2c41..b62784707 100644
--- a/numpy/core/src/multiarray/datetime.c
+++ b/numpy/core/src/multiarray/datetime.c
@@ -575,6 +575,12 @@ convert_datetime_to_datetimestruct(PyArray_DatetimeMetaData *meta,
     out->year = 1970;
     out->month = 1;
     out->day = 1;
+
+    /* NaT is signaled in the year */
+    if (dt == NPY_DATETIME_NAT) {
+        out->year = NPY_DATETIME_NAT;
+        return 0;
+    }
     
     /* Extract the event number */
     if (meta->events > 1) {
@@ -1541,6 +1547,108 @@ get_datetime_conversion_factor(PyArray_DatetimeMetaData *src_meta,
     *out_denom = (npy_int64)(denom / gcd);
 }
 
+/*
+ * Determines whether the 'divisor' metadata divides evenly into
+ * the 'dividend' metadata.
+ */
+NPY_NO_EXPORT npy_bool
+datetime_metadata_divides(
+                        PyArray_Descr *dividend,
+                        PyArray_Descr *divisor,
+                        int strict_with_nonlinear_units)
+{
+    PyArray_DatetimeMetaData *meta1, *meta2;
+    npy_uint64 num1, num2;
+
+    /* Must be datetime types */
+    if ((dividend->type_num != NPY_DATETIME &&
+                        dividend->type_num != NPY_TIMEDELTA) ||
+                    (divisor->type_num != NPY_DATETIME &&
+                        divisor->type_num != NPY_TIMEDELTA)) {
+        return 0;
+    }
+
+    meta1 = get_datetime_metadata_from_dtype(dividend);
+    if (meta1 == NULL) {
+        PyErr_Clear();
+        return 0;
+    }
+    meta2 = get_datetime_metadata_from_dtype(divisor);
+    if (meta2 == NULL) {
+        PyErr_Clear();
+        return 0;
+    }
+
+    /* Events must match */
+    if (meta1->events != meta2->events) {
+        return 0;
+    }
+
+    num1 = (npy_uint64)meta1->num;
+    num2 = (npy_uint64)meta2->num;
+
+    /* If the bases are different, factor in a conversion */
+    if (meta1->base != meta2->base) {
+        /*
+         * Years, Months, and Business days are incompatible with
+         * all other units (except years and months are compatible
+         * with each other).
+         */
+        if (meta1->base == NPY_FR_B || meta2->base == NPY_FR_B) {
+            return 0;
+        }
+        else if (meta1->base == NPY_FR_Y) {
+            if (meta2->base == NPY_FR_M) {
+                num1 *= 12;
+            }
+            else if (strict_with_nonlinear_units) {
+                return 0;
+            }
+            else {
+                /* Could do something complicated here */
+                return 1;
+            }
+        }
+        else if (meta2->base == NPY_FR_Y) {
+            if (meta1->base == NPY_FR_M) {
+                num2 *= 12;
+            }
+            else if (strict_with_nonlinear_units) {
+                return 0;
+            }
+            else {
+                /* Could do something complicated here */
+                return 1;
+            }
+        }
+        else if (meta1->base == NPY_FR_M || meta2->base == NPY_FR_M) {
+            if (strict_with_nonlinear_units) {
+                return 0;
+            }
+            else {
+                /* Could do something complicated here */
+                return 1;
+            }
+        }
+
+        /* Take the greater base (unit sizes are decreasing in enum) */
+        if (meta1->base > meta2->base) {
+            num2 *= get_datetime_units_factor(meta2->base, meta1->base);
+            if (num2 == 0) {
+                return 0;
+            }
+        }
+        else {
+            num1 *= get_datetime_units_factor(meta1->base, meta2->base);
+            if (num1 == 0) {
+                return 0;
+            }
+        }
+    }
+
+    return (num2 % num1) == 0;
+}
+
 
 NPY_NO_EXPORT PyObject *
 compute_datetime_metadata_greatest_common_divisor(
@@ -1555,8 +1663,8 @@ compute_datetime_metadata_greatest_common_divisor(
 
     if ((type1->type_num != NPY_DATETIME &&
                         type1->type_num != NPY_TIMEDELTA) ||
-                    (type1->type_num != NPY_DATETIME &&
-                        type1->type_num != NPY_TIMEDELTA)) {
+                    (type2->type_num != NPY_DATETIME &&
+                        type2->type_num != NPY_TIMEDELTA)) {
         PyErr_SetString(PyExc_TypeError,
                 "Require datetime types for metadata "
                 "greatest common divisor operation");
diff --git a/numpy/core/src/multiarray/scalartypes.c.src b/numpy/core/src/multiarray/scalartypes.c.src
index ebcfd4f96..d6ee85df0 100644
--- a/numpy/core/src/multiarray/scalartypes.c.src
+++ b/numpy/core/src/multiarray/scalartypes.c.src
@@ -3426,18 +3426,21 @@ initialize_casting_tables(void)
     for (i = 0; i < NPY_NTYPES; ++i) {
         /* Identity */
         _npy_can_cast_safely_table[i][i] = 1;
-        /* Bool -> <Anything> */
-        _npy_can_cast_safely_table[NPY_BOOL][i] = 1;
-        /* DateTime sits out for these... */
-        if (i != PyArray_DATETIME && i != PyArray_TIMEDELTA) {
-            /* <Anything> -> Object */
-            _npy_can_cast_safely_table[i][NPY_OBJECT] = 1;
-            /* <Anything> -> Void */
-            _npy_can_cast_safely_table[i][NPY_VOID] = 1;
+        if (i != NPY_DATETIME) {
+            /*
+             * Bool -> <Anything> except datetime (since
+             *                    it conceptually has no zero)
+             */
+            _npy_can_cast_safely_table[NPY_BOOL][i] = 1;
         }
+        /* <Anything> -> Object */
+        _npy_can_cast_safely_table[i][NPY_OBJECT] = 1;
+        /* <Anything> -> Void */
+        _npy_can_cast_safely_table[i][NPY_VOID] = 1;
     }
 
     _npy_can_cast_safely_table[NPY_STRING][NPY_UNICODE] = 1;
+    _npy_can_cast_safely_table[NPY_BOOL][NPY_TIMEDELTA] = 1;
 
 #ifndef NPY_SIZEOF_BYTE
 #define NPY_SIZEOF_BYTE 1
@@ -3467,8 +3470,13 @@ initialize_casting_tables(void)
 #define _FROM_BSIZE NPY_SIZEOF_@FROM_BASENAME@
 #define _FROM_NUM   (NPY_@FROM_NAME@)
 
-    _npy_can_cast_safely_table[_FROM_NUM][PyArray_STRING] = 1;
-    _npy_can_cast_safely_table[_FROM_NUM][PyArray_UNICODE] = 1;
+    _npy_can_cast_safely_table[_FROM_NUM][NPY_STRING] = 1;
+    _npy_can_cast_safely_table[_FROM_NUM][NPY_UNICODE] = 1;
+
+    /* Allow casts from any integer to the TIMEDELTA type */
+#if @from_isint@ || @from_isuint@
+    _npy_can_cast_safely_table[_FROM_NUM][NPY_TIMEDELTA] = 1;
+#endif
 
 /**begin repeat1
  * #TO_NAME = BYTE, UBYTE, SHORT, USHORT, INT, UINT, LONG, ULONG,
@@ -3651,9 +3659,6 @@ initialize_casting_tables(void)
             }
         }
     }
-    /* Special case date-time */
-    _npy_type_promotion_table[NPY_DATETIME][NPY_TIMEDELTA] = NPY_DATETIME;
-    _npy_type_promotion_table[NPY_TIMEDELTA][NPY_DATETIME] = NPY_DATETIME;
 }
 
 
diff --git a/numpy/core/src/umath/loops.c.src b/numpy/core/src/umath/loops.c.src
index 60383bc45..294eb7f47 100644
--- a/numpy/core/src/umath/loops.c.src
+++ b/numpy/core/src/umath/loops.c.src
@@ -946,7 +946,12 @@ TIMEDELTA_negative(char **args, intp *dimensions, intp *steps, void *NPY_UNUSED(
 {
     UNARY_LOOP {
         const npy_timedelta in1 = *(npy_timedelta *)ip1;
-        *((npy_timedelta *)op1) = (npy_timedelta)(-(npy_timedelta)in1);
+        if (in1 == NPY_DATETIME_NAT) {
+            *((npy_timedelta *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_timedelta *)op1) = -in1;
+        }
     }
 }
 
@@ -955,7 +960,12 @@ TIMEDELTA_absolute(char **args, intp *dimensions, intp *steps, void *NPY_UNUSED(
 {
     UNARY_LOOP {
         const npy_timedelta in1 = *(npy_timedelta *)ip1;
-        *((npy_timedelta *)op1) = (in1 >= 0) ? in1 : -in1;
+        if (in1 == NPY_DATETIME_NAT) {
+            *((npy_timedelta *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_timedelta *)op1) = (in1 >= 0) ? in1 : -in1;
+        }
     }
 }
 
@@ -1027,7 +1037,7 @@ NPY_NO_EXPORT void
     BINARY_LOOP {
         const @type@ in1 = *(@type@ *)ip1;
         const @type@ in2 = *(@type@ *)ip2;
-        *((Bool *)op1) = in1 @OP@ in2;
+        *((npy_bool *)op1) = in1 @OP@ in2;
     }
 }
 /**end repeat1**/
@@ -1042,16 +1052,24 @@ NPY_NO_EXPORT void
     if (IS_BINARY_REDUCE) {
         BINARY_REDUCE_LOOP(@type@) {
             const @type@ in2 = *(@type@ *)ip2;
-            io1 = (io1 @OP@ in2) ? io1 : in2;
+            io1 = (io1 @OP@ in2 || in2 == NPY_DATETIME_NAT) ? io1 : in2;
         }
         *((@type@ *)iop1) = io1;
     }
     else {
-    BINARY_LOOP {
-        const @type@ in1 = *(@type@ *)ip1;
-        const @type@ in2 = *(@type@ *)ip2;
-        *((@type@ *)op1) = (in1 @OP@ in2) ? in1 : in2;
-    }
+        BINARY_LOOP {
+            const @type@ in1 = *(@type@ *)ip1;
+            const @type@ in2 = *(@type@ *)ip2;
+            if (in1 == NPY_DATETIME_NAT) {
+                *((@type@ *)op1) = in2;
+            }
+            else if (in2 == NPY_DATETIME_NAT) {
+                *((@type@ *)op1) = in1;
+            }
+            else {
+                *((@type@ *)op1) = (in1 @OP@ in2) ? in1 : in2;
+            }
+        }
     }
 }
 /**end repeat1**/
@@ -1064,7 +1082,12 @@ DATETIME_Mm_M_add(char **args, intp *dimensions, intp *steps, void *NPY_UNUSED(d
     BINARY_LOOP {
         const datetime in1 = *(datetime *)ip1;
         const timedelta in2 = *(timedelta *)ip2;
-        *((datetime *)op1) = in1 + in2;
+        if (in1 == NPY_DATETIME_NAT || in2 == NPY_DATETIME_NAT) {
+            *((npy_datetime *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_datetime *)op1) = in1 + in2;
+        }
     }
 }
 
@@ -1074,7 +1097,12 @@ DATETIME_mM_M_add(char **args, intp *dimensions, intp *steps, void *NPY_UNUSED(f
     BINARY_LOOP {
         const timedelta in1 = *(timedelta *)ip1;
         const datetime in2 = *(datetime *)ip2;
-        *((datetime *)op1) = in1 + in2;
+        if (in1 == NPY_DATETIME_NAT || in2 == NPY_DATETIME_NAT) {
+            *((npy_datetime *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_datetime *)op1) = in1 + in2;
+        }
     }
 }
 
@@ -1084,7 +1112,12 @@ TIMEDELTA_mm_m_add(char **args, intp *dimensions, intp *steps, void *NPY_UNUSED(
     BINARY_LOOP {
         const timedelta in1 = *(timedelta *)ip1;
         const timedelta in2 = *(timedelta *)ip2;
-        *((timedelta *)op1) = in1 + in2;
+        if (in1 == NPY_DATETIME_NAT || in2 == NPY_DATETIME_NAT) {
+            *((npy_timedelta *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_timedelta *)op1) = in1 + in2;
+        }
     }
 }
 
@@ -1094,7 +1127,12 @@ DATETIME_Mm_M_subtract(char **args, intp *dimensions, intp *steps, void *NPY_UNU
     BINARY_LOOP {
         const datetime in1 = *(datetime *)ip1;
         const timedelta in2 = *(timedelta *)ip2;
-        *((datetime *)op1) = in1 - in2;
+        if (in1 == NPY_DATETIME_NAT || in2 == NPY_DATETIME_NAT) {
+            *((npy_datetime *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_datetime *)op1) = in1 - in2;
+        }
     }
 }
 
@@ -1104,7 +1142,12 @@ DATETIME_MM_m_subtract(char **args, intp *dimensions, intp *steps, void *NPY_UNU
     BINARY_LOOP {
         const datetime in1 = *(datetime *)ip1;
         const datetime in2 = *(datetime *)ip2;
-        *((timedelta *)op1) = in1 - in2;
+        if (in1 == NPY_DATETIME_NAT || in2 == NPY_DATETIME_NAT) {
+            *((npy_timedelta *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_timedelta *)op1) = in1 - in2;
+        }
     }
 }
 
@@ -1114,7 +1157,12 @@ TIMEDELTA_mm_m_subtract(char **args, intp *dimensions, intp *steps, void *NPY_UN
     BINARY_LOOP {
         const timedelta in1 = *(timedelta *)ip1;
         const timedelta in2 = *(timedelta *)ip2;
-        *((timedelta *)op1) = in1 - in2;
+        if (in1 == NPY_DATETIME_NAT || in2 == NPY_DATETIME_NAT) {
+            *((npy_timedelta *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_timedelta *)op1) = in1 - in2;
+        }
     }
 }
 
@@ -1125,7 +1173,12 @@ TIMEDELTA_mq_m_multiply(char **args, intp *dimensions, intp *steps, void *NPY_UN
     BINARY_LOOP {
         const npy_timedelta in1 = *(npy_timedelta *)ip1;
         const npy_int64 in2 = *(npy_int64 *)ip2;
-        *((npy_timedelta *)op1) = in1 * in2;
+        if (in1 == NPY_DATETIME_NAT) {
+            *((npy_timedelta *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_timedelta *)op1) = in1 * in2;
+        }
     }
 }
 
@@ -1136,7 +1189,12 @@ TIMEDELTA_qm_m_multiply(char **args, intp *dimensions, intp *steps, void *NPY_UN
     BINARY_LOOP {
         const npy_int64 in1 = *(npy_int64 *)ip1;
         const npy_timedelta in2 = *(npy_timedelta *)ip2;
-        *((npy_timedelta *)op1) = in1 * in2;
+        if (in2 == NPY_DATETIME_NAT) {
+            *((npy_timedelta *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_timedelta *)op1) = in1 * in2;
+        }
     }
 }
 
@@ -1146,7 +1204,12 @@ TIMEDELTA_md_m_multiply(char **args, intp *dimensions, intp *steps, void *NPY_UN
     BINARY_LOOP {
         const npy_timedelta in1 = *(npy_timedelta *)ip1;
         const double in2 = *(double *)ip2;
-        *((npy_timedelta *)op1) = (npy_timedelta)(in1 * in2);
+        if (in1 == NPY_DATETIME_NAT || isnan(in2)) {
+            *((npy_timedelta *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_timedelta *)op1) = (npy_timedelta)(in1 * in2);
+        }
     }
 }
 
@@ -1156,7 +1219,12 @@ TIMEDELTA_dm_m_multiply(char **args, intp *dimensions, intp *steps, void *NPY_UN
     BINARY_LOOP {
         const double in1 = *(double *)ip1;
         const npy_timedelta in2 = *(npy_timedelta *)ip2;
-        *((npy_timedelta *)op1) = (npy_timedelta)(in1 * in2);
+        if (isnan(in1) || in2 == NPY_DATETIME_NAT) {
+            *((npy_timedelta *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_timedelta *)op1) = (npy_timedelta)(in1 * in2);
+        }
     }
 }
 
@@ -1167,7 +1235,12 @@ TIMEDELTA_mq_m_divide(char **args, intp *dimensions, intp *steps, void *NPY_UNUS
     BINARY_LOOP {
         const npy_timedelta in1 = *(npy_timedelta *)ip1;
         const npy_int64 in2 = *(npy_int64 *)ip2;
-        *((npy_timedelta *)op1) = in1 / in2;
+        if (in1 == NPY_DATETIME_NAT || in2 == 0) {
+            *((npy_timedelta *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            *((npy_timedelta *)op1) = in1 / in2;
+        }
     }
 }
 
@@ -1177,7 +1250,18 @@ TIMEDELTA_md_m_divide(char **args, intp *dimensions, intp *steps, void *NPY_UNUS
     BINARY_LOOP {
         const npy_timedelta in1 = *(npy_timedelta *)ip1;
         const double in2 = *(double *)ip2;
-        *((npy_timedelta *)op1) = (npy_timedelta)(in1 / in2);
+        if (in1 == NPY_DATETIME_NAT) {
+            *((npy_timedelta *)op1) = NPY_DATETIME_NAT;
+        }
+        else {
+            double result = in1 / in2;
+            if (isnan(result)) {
+                *((npy_timedelta *)op1) = NPY_DATETIME_NAT;
+            }
+            else {
+                *((npy_timedelta *)op1) = (npy_timedelta)(result);
+            }
+        }
     }
 }
 
diff --git a/numpy/core/src/umath/ufunc_object.c b/numpy/core/src/umath/ufunc_object.c
index 95c8463ad..b9e4c302e 100644
--- a/numpy/core/src/umath/ufunc_object.c
+++ b/numpy/core/src/umath/ufunc_object.c
@@ -2164,7 +2164,8 @@ PyUFunc_AdditionTypeResolution(PyUFuncObject *ufunc,
             Py_INCREF(out_dtypes[2]);
         }
         /* m8[<A>] + int => m8[<A>] + m8[<A>] */
-        else if (PyTypeNum_ISINTEGER(type_num2)) {
+        else if (PyTypeNum_ISINTEGER(type_num2) ||
+                                    PyTypeNum_ISBOOL(type_num2)) {
             out_dtypes[0] = PyArray_DESCR(operands[0]);
             Py_INCREF(out_dtypes[0]);
             out_dtypes[1] = out_dtypes[0];
@@ -2182,7 +2183,8 @@ PyUFunc_AdditionTypeResolution(PyUFuncObject *ufunc,
         /* M8[<A>] + m8[<B>] => M8[<A>] + m8[<A>] */
         /* M8[<A>] + int => M8[<A>] + m8[<A>] */
         if (type_num2 == NPY_TIMEDELTA ||
-                    PyTypeNum_ISINTEGER(type_num2)) {
+                    PyTypeNum_ISINTEGER(type_num2) ||
+                    PyTypeNum_ISBOOL(type_num2)) {
             /* Make a new NPY_TIMEDELTA, and copy type1's metadata */
             out_dtypes[1] = timedelta_dtype_with_copied_meta(
                                             PyArray_DESCR(operands[0]));
@@ -2200,7 +2202,7 @@ PyUFunc_AdditionTypeResolution(PyUFuncObject *ufunc,
             goto type_reso_error;
         }
     }
-    else if (PyTypeNum_ISINTEGER(type_num1)) {
+    else if (PyTypeNum_ISINTEGER(type_num1) || PyTypeNum_ISBOOL(type_num1)) {
         /* int + m8[<A>] => m8[<A>] + m8[<A>] */
         if (type_num2 == NPY_TIMEDELTA) {
             out_dtypes[0] = PyArray_DESCR(operands[1]);
@@ -2325,7 +2327,8 @@ PyUFunc_SubtractionTypeResolution(PyUFuncObject *ufunc,
             Py_INCREF(out_dtypes[2]);
         }
         /* m8[<A>] - int => m8[<A>] - m8[<A>] */
-        else if (PyTypeNum_ISINTEGER(type_num2)) {
+        else if (PyTypeNum_ISINTEGER(type_num2) ||
+                                        PyTypeNum_ISBOOL(type_num2)) {
             out_dtypes[0] = PyArray_DESCR(operands[0]);
             Py_INCREF(out_dtypes[0]);
             out_dtypes[1] = out_dtypes[0];
@@ -2343,7 +2346,8 @@ PyUFunc_SubtractionTypeResolution(PyUFuncObject *ufunc,
         /* M8[<A>] - m8[<B>] => M8[<A>] - m8[<A>] */
         /* M8[<A>] - int => M8[<A>] - m8[<A>] */
         if (type_num2 == NPY_TIMEDELTA ||
-                    PyTypeNum_ISINTEGER(type_num2)) {
+                    PyTypeNum_ISINTEGER(type_num2) ||
+                    PyTypeNum_ISBOOL(type_num2)) {
             /* Make a new NPY_TIMEDELTA, and copy type1's metadata */
             out_dtypes[1] = timedelta_dtype_with_copied_meta(
                                             PyArray_DESCR(operands[0]));
@@ -2395,7 +2399,7 @@ PyUFunc_SubtractionTypeResolution(PyUFuncObject *ufunc,
             goto type_reso_error;
         }
     }
-    else if (PyTypeNum_ISINTEGER(type_num1)) {
+    else if (PyTypeNum_ISINTEGER(type_num1) || PyTypeNum_ISBOOL(type_num1)) {
         /* int - m8[<A>] => m8[<A>] - m8[<A>] */
         if (type_num2 == NPY_TIMEDELTA) {
             out_dtypes[0] = PyArray_DESCR(operands[1]);
@@ -2491,7 +2495,7 @@ PyUFunc_MultiplicationTypeResolution(PyUFuncObject *ufunc,
 
     if (type_num1 == NPY_TIMEDELTA) {
         /* m8[<A>] * int## => m8[<A>] * int64 */
-        if (PyTypeNum_ISINTEGER(type_num2)) {
+        if (PyTypeNum_ISINTEGER(type_num2) || PyTypeNum_ISBOOL(type_num2)) {
             out_dtypes[0] = PyArray_DESCR(operands[0]);
             Py_INCREF(out_dtypes[0]);
             out_dtypes[1] = PyArray_DescrNewFromType(NPY_INT64);
@@ -2524,7 +2528,7 @@ PyUFunc_MultiplicationTypeResolution(PyUFuncObject *ufunc,
             goto type_reso_error;
         }
     }
-    else if (PyTypeNum_ISINTEGER(type_num1)) {
+    else if (PyTypeNum_ISINTEGER(type_num1) || PyTypeNum_ISBOOL(type_num1)) {
         /* int## * m8[<A>] => int64 * m8[<A>] */
         if (type_num2 == NPY_TIMEDELTA) {
             out_dtypes[0] = PyArray_DescrNewFromType(NPY_INT64);
diff --git a/numpy/core/tests/test_datetime.py b/numpy/core/tests/test_datetime.py
index 5b8aff14f..82681d58f 100644
--- a/numpy/core/tests/test_datetime.py
+++ b/numpy/core/tests/test_datetime.py
@@ -5,7 +5,7 @@ from numpy.compat import asbytes
 import datetime
 
 class TestDateTime(TestCase):
-    def test_creation(self):
+    def test_datetime_dtype_creation(self):
         for unit in ['Y', 'M', 'W', 'B', 'D',
                      'h', 'm', 's', 'ms', 'us',
                      'ns', 'ps', 'fs', 'as']:
@@ -26,6 +26,60 @@ class TestDateTime(TestCase):
         assert_raises(TypeError, np.dtype, 'M16')
         assert_raises(TypeError, np.dtype, 'm16')
 
+    def test_datetime_casting_rules(self):
+        # Cannot cast safely/same_kind between timedelta and datetime
+        assert_(not np.can_cast('m8', 'M8', casting='same_kind'))
+        assert_(not np.can_cast('M8', 'm8', casting='same_kind'))
+        assert_(not np.can_cast('m8', 'M8', casting='safe'))
+        assert_(not np.can_cast('M8', 'm8', casting='safe'))
+
+        # Can cast safely/same_kind from integer to timedelta
+        assert_(np.can_cast('i8', 'm8', casting='same_kind'))
+        assert_(np.can_cast('i8', 'm8', casting='safe'))
+
+        # Cannot cast safely/same_kind from float to timedelta
+        assert_(not np.can_cast('f4', 'm8', casting='same_kind'))
+        assert_(not np.can_cast('f4', 'm8', casting='safe'))
+
+        # Cannot cast safely/same_kind from integer to datetime
+        assert_(not np.can_cast('i8', 'M8', casting='same_kind'))
+        assert_(not np.can_cast('i8', 'M8', casting='safe'))
+
+        # Cannot cast safely/same_kind from bool to datetime
+        assert_(not np.can_cast('b1', 'M8', casting='same_kind'))
+        assert_(not np.can_cast('b1', 'M8', casting='safe'))
+        # Can cast safely/same_kind from bool to timedelta
+        assert_(np.can_cast('b1', 'm8', casting='same_kind'))
+        assert_(np.can_cast('b1', 'm8', casting='safe'))
+
+        # Can cast datetime safely from months/years to days
+        assert_(np.can_cast('M8[M]', 'M8[D]', casting='safe'))
+        assert_(np.can_cast('M8[Y]', 'M8[D]', casting='safe'))
+        # Cannot cast timedelta safely from months/years to days
+        assert_(not np.can_cast('m8[M]', 'm8[D]', casting='safe'))
+        assert_(not np.can_cast('m8[Y]', 'm8[D]', casting='safe'))
+        # Can cast same_kind from months/years to days
+        assert_(np.can_cast('M8[M]', 'M8[D]', casting='same_kind'))
+        assert_(np.can_cast('M8[Y]', 'M8[D]', casting='same_kind'))
+        assert_(np.can_cast('m8[M]', 'm8[D]', casting='same_kind'))
+        assert_(np.can_cast('m8[Y]', 'm8[D]', casting='same_kind'))
+
+        # Cannot cast safely if the integer multiplier doesn't divide
+        assert_(not np.can_cast('M8[7h]', 'M8[3h]', casting='safe'))
+        assert_(not np.can_cast('M8[3h]', 'M8[6h]', casting='safe'))
+        # But can cast same_kind
+        assert_(np.can_cast('M8[7h]', 'M8[3h]', casting='same_kind'))
+        # Can cast safely if the integer multiplier does divide
+        assert_(np.can_cast('M8[6h]', 'M8[3h]', casting='safe'))
+
+    def test_datetime_nat_casting(self):
+        a = np.array('NaT', dtype='M8[D]')
+        assert_equal(a.astype('M8[s]'), np.array('NaT', dtype='M8[s]'))
+        assert_equal(a.astype('M8[ms]'), np.array('NaT', dtype='M8[ms]'))
+        assert_equal(a.astype('M8[M]'), np.array('NaT', dtype='M8[M]'))
+        assert_equal(a.astype('M8[Y]'), np.array('NaT', dtype='M8[Y]'))
+        assert_equal(a.astype('M8[W]'), np.array('NaT', dtype='M8[W]'))
+
     def test_days_creation(self):
         assert_equal(np.array('1599', dtype='M8[D]').astype('i8'),
                 (1600-1970)*365 - (1972-1600)/4 + 3 - 365)
@@ -163,7 +217,7 @@ class TestDateTime(TestCase):
 
     def test_pyobject_roundtrip(self):
         # All datetime types should be able to roundtrip through object
-        a = np.array([0,0,0,0,0,0,0,0,
+        a = np.array([0,0,0,0,0,0,0,0,0,
                       -1020040340, -2942398, -1, 0, 1, 234523453, 1199164176],
                                                         dtype=np.int64)
         for unit in ['M8[as]', 'M8[16fs]', 'M8[ps]', 'M8[us]',
@@ -178,6 +232,7 @@ class TestDateTime(TestCase):
             b[5] = '1970-01-01'
             b[6] = '9999-12-31T23:59:59.999999Z'
             b[7] = '10000-01-01'
+            b[8] = 'NaT'
 
             assert_equal(b.astype(object).astype(unit), b,
                             "Error roundtripping unit %s" % unit)
@@ -222,10 +277,21 @@ class TestDateTime(TestCase):
                 assert_equal(np.array('today', dtype=dt1),
                              np.array('today', dtype=dt2))
 
+        # Shouldn't be able to compare datetime and timedelta
+        # TODO: Changing to 'same_kind' or 'safe' casting in the ufuncs by
+        #       default is needed to properly catch this kind of thing...
+        #a = np.array('2012-12-21', dtype='M8[D]')
+        #b = np.array(3, dtype='m8[D]')
+        #assert_raises(TypeError, np.less, a, b)
+        #a = np.array('2012-12-21', dtype='M8[D]')
+        #b = np.array(3, dtype='m8[D]')
+        #assert_raises(TypeError, np.less, a, b, casting='same_kind')
+
     def test_datetime_add(self):
         dta = np.array('2012-12-21', dtype='M8[D]')
         dtb = np.array('2012-12-24', dtype='M8[D]')
         dtc = np.array('1940-12-24', dtype='M8[D]')
+        dtnat = np.array('NaT', dtype='M8[D]')
         tda = np.array(3, dtype='m8[D]')
         tdb = np.array(11, dtype='m8[h]')
         tdc = np.array(3*24 + 11, dtype='m8[h]')
@@ -233,23 +299,41 @@ class TestDateTime(TestCase):
         # m8 + m8
         assert_equal(tda + tdb, tdc)
         assert_equal((tda + tdb).dtype, np.dtype('m8[h]'))
+        # m8 + bool
+        assert_equal(tdb + True, tdb + 1)
+        assert_equal((tdb + True).dtype, np.dtype('m8[h]'))
         # m8 + int
         assert_equal(tdb + 3*24, tdc)
         assert_equal((tdb + 3*24).dtype, np.dtype('m8[h]'))
+        # bool + m8
+        assert_equal(False + tdb, tdb)
+        assert_equal((False + tdb).dtype, np.dtype('m8[h]'))
         # int + m8
         assert_equal(3*24 + tdb, tdc)
         assert_equal((3*24 + tdb).dtype, np.dtype('m8[h]'))
+        # M8 + bool
+        assert_equal(dta + True, dta + 1)
+        assert_equal(dtnat + True, dtnat)
+        assert_equal((dta + True).dtype, np.dtype('M8[D]'))
         # M8 + int
         assert_equal(dta + 3, dtb)
+        assert_equal(dtnat + 3, dtnat)
         assert_equal((dta + 3).dtype, np.dtype('M8[D]'))
+        # bool + M8
+        assert_equal(False + dta, dta)
+        assert_equal(False + dtnat, dtnat)
+        assert_equal((False + dta).dtype, np.dtype('M8[D]'))
         # int + M8
         assert_equal(3 + dta, dtb)
+        assert_equal(3 + dtnat, dtnat)
         assert_equal((3 + dta).dtype, np.dtype('M8[D]'))
         # M8 + m8
         assert_equal(dta + tda, dtb)
+        assert_equal(dtnat + tda, dtnat)
         assert_equal((dta + tda).dtype, np.dtype('M8[D]'))
         # m8 + M8
         assert_equal(tda + dta, dtb)
+        assert_equal(tda + dtnat, dtnat)
         assert_equal((tda + dta).dtype, np.dtype('M8[D]'))
 
         # In M8 + m8, the M8 controls the result type
@@ -270,6 +354,7 @@ class TestDateTime(TestCase):
         dtb = np.array('2012-12-24', dtype='M8[D]')
         dtc = np.array('1940-12-24', dtype='M8[D]')
         dtd = np.array('1940-12-24', dtype='M8[h]')
+        dtnat = np.array('NaT', dtype='M8[D]')
         tda = np.array(3, dtype='m8[D]')
         tdb = np.array(11, dtype='m8[h]')
         tdc = np.array(3*24 - 11, dtype='m8[h]')
@@ -279,17 +364,29 @@ class TestDateTime(TestCase):
         assert_equal((tda - tdb).dtype, np.dtype('m8[h]'))
         assert_equal(tdb - tda, -tdc)
         assert_equal((tdb - tda).dtype, np.dtype('m8[h]'))
+        # m8 - bool
+        assert_equal(tdc - True, tdc - 1)
+        assert_equal((tdc - True).dtype, np.dtype('m8[h]'))
         # m8 - int
         assert_equal(tdc - 3*24, -tdb)
         assert_equal((tdc - 3*24).dtype, np.dtype('m8[h]'))
         # int - m8
+        assert_equal(False - tdb, -tdb)
+        assert_equal((False - tdb).dtype, np.dtype('m8[h]'))
+        # int - m8
         assert_equal(3*24 - tdb, tdc)
         assert_equal((3*24 - tdb).dtype, np.dtype('m8[h]'))
+        # M8 - bool
+        assert_equal(dtb - True, dtb - 1)
+        assert_equal(dtnat - True, dtnat)
+        assert_equal((dtb - True).dtype, np.dtype('M8[D]'))
         # M8 - int
         assert_equal(dtb - 3, dta)
+        assert_equal(dtnat - 3, dtnat)
         assert_equal((dtb - 3).dtype, np.dtype('M8[D]'))
         # M8 - m8
         assert_equal(dtb - tda, dta)
+        assert_equal(dtnat - tda, dtnat)
         assert_equal((dtb - tda).dtype, np.dtype('M8[D]'))
 
         # In M8 - m8, the M8 controls the result type
@@ -302,6 +399,8 @@ class TestDateTime(TestCase):
         assert_raises(TypeError, np.subtract, dtc, dtd)
         # m8 - M8
         assert_raises(TypeError, np.subtract, tda, dta)
+        # bool - M8
+        assert_raises(TypeError, np.subtract, False, dta)
         # int - M8
         assert_raises(TypeError, np.subtract, 3, dta)
 
@@ -388,6 +487,14 @@ class TestDateTime(TestCase):
         # of the units chosen
         assert_equal(np.minimum(a.view('i8'),b.view('i8')), a.view('i8'))
 
+        # Interaction with NaT
+        a = np.array('1999-03-12T13Z', dtype='M8[2m]')
+        dtnat = np.array('NaT', dtype='M8[D]')
+        assert_equal(np.minimum(a,dtnat), a)
+        assert_equal(np.minimum(dtnat,a), a)
+        assert_equal(np.maximum(a,dtnat), a)
+        assert_equal(np.maximum(dtnat,a), a)
+
         # Also do timedelta
         a = np.array(3, dtype='m8[h]')
         b = np.array(3*3600 - 3, dtype='m8[s]')
@@ -403,6 +510,17 @@ class TestDateTime(TestCase):
         # of the units chosen
         assert_equal(np.minimum(a.view('i8'),b.view('i8')), a.view('i8'))
 
+        # should raise between datetime and timedelta
+        #
+        # TODO: Allowing unsafe casting by
+        #       default in ufuncs strikes again... :(
+        #a = np.array(3, dtype='m8[h]')
+        #b = np.array('1999-03-12T12Z', dtype='M8[s]')
+        #assert_raises(TypeError, np.minimum, a, b)
+        #assert_raises(TypeError, np.maximum, a, b)
+        #assert_raises(TypeError, np.fmin, a, b)
+        #assert_raises(TypeError, np.fmax, a, b)
+
     def test_hours(self):
         t = np.ones(3, dtype='M8[s]')
         t[0] = 60*60*24 + 60*60*10