Merge tree-ssa-20020619-branch into mainline.

git-svn-id: svn+ssh://gcc.gnu.org/svn/gcc/trunk@81764 138bc75d-0d04-0410-961f-82ee72b054a4

18 files changed, 2520 insertions, 0 deletions

diff --git a/libgfortran/intrinsics/abort.c b/libgfortran/intrinsics/abort.c
new file mode 100644
index 00000000000..26a25362a7d
--- /dev/null
+++ b/libgfortran/intrinsics/abort.c
@@ -0,0 +1,31 @@
+/* Implementation of the ABORT intrinsic.
+   Copyright (C) 2003 Free Software Foundation, Inc.
+
+This file is part of the GNU Fortran 95 runtime library (libgfortran).
+
+Libgfor is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Libgfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "libgfortran.h"
+#include <assert.h>
+
+
+void prefix(abort) (void);
+
+void prefix(abort) ()
+{
+  abort ();
+}
+
diff --git a/libgfortran/intrinsics/associated.c b/libgfortran/intrinsics/associated.c
new file mode 100644
index 00000000000..ba52a205d12
--- /dev/null
+++ b/libgfortran/intrinsics/associated.c
@@ -0,0 +1,50 @@
+/* Implementation of the ASSOCIATED intrinsic
+   Copyright 2003 Free Software Foundation, Inc.
+   Contributed by kejia Zhao (CCRG) <kejia_zh@yahoo.com.cn>
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfor is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Ligbfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "libgfortran.h"
+
+#define associated prefix(associated)
+
+enum { FALSE = 0, TRUE = 1 };
+
+
+GFC_LOGICAL_4
+associated (const gfc_array_void *pointer, const gfc_array_void *target)
+{
+  int n, rank;
+
+  if (GFC_DESCRIPTOR_DATA (pointer) != GFC_DESCRIPTOR_DATA (target))
+    return FALSE;
+  if (GFC_DESCRIPTOR_DTYPE (pointer) != GFC_DESCRIPTOR_DTYPE (target))
+    return FALSE;
+
+  rank = GFC_DESCRIPTOR_RANK (pointer);
+  for (n = 0; n < rank; n++)
+    {
+      if (pointer->dim[n].stride != target->dim[n].stride)
+        return FALSE;
+      if ((pointer->dim[n].ubound - pointer->dim[n].lbound)
+          != (target->dim[n].ubound - target->dim[n].lbound))
+        return FALSE;
+    }
+
+  return  TRUE;
+}
diff --git a/libgfortran/intrinsics/cpu_time.c b/libgfortran/intrinsics/cpu_time.c
new file mode 100644
index 00000000000..9fd954b9a71
--- /dev/null
+++ b/libgfortran/intrinsics/cpu_time.c
@@ -0,0 +1,116 @@
+/* Implementation of the CPU_TIME intrinsic.
+   Copyright (C) 2003 Free Software Foundation, Inc.
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfor is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Libgfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "config.h"
+#include <sys/types.h>
+#include "libgfortran.h"
+
+#ifdef HAVE_UNISTD_H
+#include <unistd.h>
+#endif
+
+/* The CPU_TIME intrinsic to "compare different algorithms on the same
+   computer or discover which parts are the most expensive", so we
+   need a way to get the CPU time with the finest resolution possible.
+   We can only be accurate up to microseconds.
+
+   As usual with UNIX systems, unfortunately no single way is
+   available for all systems.  */
+
+#ifdef TIME_WITH_SYS_TIME
+#  include <sys/time.h>
+#  include <time.h>
+#else
+#  if HAVE_SYS_TIME_H
+#    include <sys/time.h>
+#  else
+#    ifdef HAVE_TIME_H
+#      include <time.h>
+#    endif
+#  endif
+#endif
+
+/* The most accurate way to get the CPU time is getrusage ().
+   If we have times(), that's good enough, too.  */
+#if defined (HAVE_GETRUSAGE) && defined (HAVE_SYS_RESOURCE_H)
+#  include <sys/resource.h>
+#else
+/* For times(), we _must_ know the number of clock ticks per second.  */
+#  if defined (HAVE_TIMES) && (defined (HZ) || defined (_SC_CLK_TCK) || defined (CLK_TCK))
+#    ifdef HAVE_SYS_PARAM_H
+#      include <sys/param.h>
+#    endif
+#    include <sys/times.h>
+#    ifndef HZ
+#      if defined _SC_CLK_TCK
+#        define HZ  sysconf(_SC_CLK_TCK)
+#      else
+#        define HZ  CLK_TCK
+#      endif
+#    endif
+#  endif  /* HAVE_TIMES etc.  */
+#endif  /* HAVE_GETRUSAGE && HAVE_SYS_RESOURCE_H  */
+
+#if defined (__GNUC__) && (__GNUC__ >= 3)
+#  define ATTRIBUTE_ALWAYS_INLINE __attribute__ ((__always_inline__))
+#else
+#  define ATTRIBUTE_ALWAYS_INLINE
+#endif
+
+static inline void __cpu_time_1 (long *, long *) ATTRIBUTE_ALWAYS_INLINE;
+
+/* Helper function for the actual implementation of the CPU_TIME
+   intrnsic.  Returns a CPU time in microseconds or -1 if no CPU time
+   could be computed.  */
+static inline void
+__cpu_time_1 (long *sec, long *usec)
+{
+#if defined (HAVE_GETRUSAGE) && defined (HAVE_SYS_RESOURCE_H)
+  struct rusage usage;
+  getrusage (0, &usage);
+  *sec = usage.ru_utime.tv_sec + usage.ru_stime.tv_sec;
+  *usec = usage.ru_utime.tv_usec + usage.ru_stime.tv_usec;
+#else /* ! HAVE_GETRUSAGE || ! HAVE_SYS_RESOURCE_H  */
+#ifdef HAVE_TIMES
+  struct tms buf;
+  times (&buf);
+  *sec = 0;
+  *usec = (buf.tms_utime + buf.tms_stime) * (1000000 / HZ);
+#else /* ! HAVE_TIMES */
+  /* We have nothing to go on.  Return -1.  */
+  *sec = -1;
+  *usec = 0;
+#endif  /* HAVE_TIMES */
+#endif  /* HAVE_GETRUSAGE */
+}
+
+#undef CPU_TIME
+#define CPU_TIME(KIND)						\
+void prefix(cpu_time_##KIND) (GFC_REAL_##KIND *__time)		\
+{								\
+  long sec, usec;						\
+  __cpu_time_1 (&sec, &usec);					\
+  *__time = (GFC_REAL_##KIND) sec +				\
+		((GFC_REAL_##KIND) usec) * 1.e-6;		\
+}
+
+CPU_TIME(4)
+CPU_TIME(8)
+
diff --git a/libgfortran/intrinsics/cshift0.c b/libgfortran/intrinsics/cshift0.c
new file mode 100644
index 00000000000..5a2c8caabe4
--- /dev/null
+++ b/libgfortran/intrinsics/cshift0.c
@@ -0,0 +1,169 @@
+/* Generic implementation of the CSHIFT intrinsic
+   Copyright 2003 Free Software Foundation, Inc.
+   Contributed by Feng Wang <wf_cs@yahoo.com>
+
+This file is part of the GNU Fortran 95 runtime library (libgfortran).
+
+Libgfortran is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Ligbfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "config.h"
+#include <stdlib.h>
+#include <assert.h>
+#include <string.h>
+#include "libgfortran.h"
+
+/* TODO: make this work for large shifts when
+   sizeof(int) < sizeof (index_type).  */
+
+static void
+__cshift0 (const gfc_array_char * ret, const gfc_array_char * array,
+    int shift, int which)
+{
+  /* r.* indicates the return array.  */
+  index_type rstride[GFC_MAX_DIMENSIONS - 1];
+  index_type rstride0;
+  index_type roffset;
+  char *rptr;
+  char *dest;
+  /* s.* indicates the source array.  */
+  index_type sstride[GFC_MAX_DIMENSIONS - 1];
+  index_type sstride0;
+  index_type soffset;
+  const char *sptr;
+  const char *src;
+
+  index_type count[GFC_MAX_DIMENSIONS - 1];
+  index_type extent[GFC_MAX_DIMENSIONS - 1];
+  index_type dim;
+  index_type size;
+  index_type len;
+  index_type n;
+
+  if (which < 1 || which > GFC_DESCRIPTOR_RANK (array))
+    runtime_error ("Argument 'DIM' is out of range in call to 'CSHIFT'");
+
+  size = GFC_DESCRIPTOR_SIZE (ret);
+
+  which = which - 1;
+
+  extent[0] = 1;
+  count[0] = 0;
+  size = GFC_DESCRIPTOR_SIZE (array);
+  n = 0;
+
+/* Initialized for avoiding compiler warnings.  */
+  roffset = size;
+  soffset = size;
+  len = 0;
+
+  for (dim = 0; dim < GFC_DESCRIPTOR_RANK (array); dim++)
+    {
+      if (dim == which)
+        {
+          roffset = ret->dim[dim].stride * size;
+          if (roffset == 0)
+            roffset = size;
+          soffset = array->dim[dim].stride * size;
+          if (soffset == 0)
+            soffset = size;
+          len = array->dim[dim].ubound + 1 - array->dim[dim].lbound;
+        }
+      else
+        {
+          count[n] = 0;
+          extent[n] = array->dim[dim].ubound + 1 - array->dim[dim].lbound;
+          rstride[n] = ret->dim[dim].stride * size;
+          sstride[n] = array->dim[dim].stride * size;
+          n++;
+        }
+    }
+  if (sstride[0] == 0)
+    sstride[0] = size;
+  if (rstride[0] == 0)
+    rstride[0] = size;
+
+  dim = GFC_DESCRIPTOR_RANK (array);
+  rstride0 = rstride[0];
+  sstride0 = sstride[0];
+  rptr = ret->data;
+  sptr = array->data;
+
+  shift = (div (shift, len)).rem;
+  if (shift < 0)
+    shift += len;
+
+  while (rptr)
+    {
+      /* Do the shift for this dimension.  */
+      src = &sptr[shift * soffset];
+      dest = rptr;
+      for (n = 0; n < len; n++)
+        {
+          memcpy (dest, src, size);
+          dest += roffset;
+          if (n == len - shift - 1)
+            src = sptr;
+          else
+            src += soffset;
+        }
+
+      /* Advance to the next section.  */
+      rptr += rstride0;
+      sptr += sstride0;
+      count[0]++;
+      n = 0;
+      while (count[n] == extent[n])
+        {
+          /* When we get to the end of a dimension, reset it and increment
+             the next dimension.  */
+          count[n] = 0;
+          /* We could precalculate these products, but this is a less
+             frequently used path so proabably not worth it.  */
+          rptr -= rstride[n] * extent[n];
+          sptr -= sstride[n] * extent[n];
+          n++;
+          if (n >= dim - 1)
+            {
+              /* Break out of the loop.  */
+              rptr = NULL;
+              break;
+            }
+          else
+            {
+              count[n]++;
+              rptr += rstride[n];
+              sptr += sstride[n];
+            }
+        }
+    }
+}
+
+
+void
+__cshift0_4 (const gfc_array_char * ret, const gfc_array_char * array,
+    const GFC_INTEGER_4 * pshift, const GFC_INTEGER_4 * pdim)
+{
+  __cshift0 (ret, array, *pshift, pdim ? *pdim : 1);
+}
+
+
+void
+__cshift0_8 (const gfc_array_char * ret, const gfc_array_char * array,
+    const GFC_INTEGER_8 * pshift, const GFC_INTEGER_8 * pdim)
+{
+  __cshift0 (ret, array, *pshift, pdim ? *pdim : 1);
+}
+
diff --git a/libgfortran/intrinsics/dprod_r8.f90 b/libgfortran/intrinsics/dprod_r8.f90
new file mode 100644
index 00000000000..d0f2063d142
--- /dev/null
+++ b/libgfortran/intrinsics/dprod_r8.f90
@@ -0,0 +1,27 @@
+!   Copyright 2003 Free Software Foundation, Inc.
+!   Contributed by Paul Brook <paul@nowt.org>
+!
+!This file is part of the GNU Fortran 95 runtime library (libgfortran).
+!
+!GNU libgfor is free software; you can redistribute it and/or
+!modify it under the terms of the GNU Lesser General Public
+!License as published by the Free Software Foundation; either
+!version 2.1 of the License, or (at your option) any later version.
+!
+!GNU libgfor is distributed in the hope that it will be useful,
+!but WITHOUT ANY WARRANTY; without even the implied warranty of
+!MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+!GNU Lesser General Public License for more details.
+!
+!You should have received a copy of the GNU Lesser General Public
+!License along with libgfor; see the file COPYING.  If not,
+!write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+!Boston, MA 02111-1307, USA.
+
+
+elemental function specific__dprod_r8 (p1, p2)
+   real (kind=4), intent (in) :: p1, p2
+   real (kind=8) :: specific__dprod_r8
+
+   specific__dprod_r8 = dprod (p1, p2)
+end function
diff --git a/libgfortran/intrinsics/eoshift0.c b/libgfortran/intrinsics/eoshift0.c
new file mode 100644
index 00000000000..f86f4bd883f
--- /dev/null
+++ b/libgfortran/intrinsics/eoshift0.c
@@ -0,0 +1,188 @@
+/* Generic implementation of the RESHAPE intrinsic
+   Copyright 2002 Free Software Foundation, Inc.
+   Contributed by Paul Brook <paul@nowt.org>
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfor is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Ligbfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "config.h"
+#include <stdlib.h>
+#include <assert.h>
+#include <string.h>
+#include "libgfortran.h"
+
+static const char zeros[16] =
+  {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
+
+/* TODO: make this work for large shifts when
+   sizeof(int) < sizeof (index_type).  */
+
+static void
+__eoshift0 (const gfc_array_char * ret, const gfc_array_char * array,
+    int shift, const char * pbound, int which)
+{
+  /* r.* indicates the return array.  */
+  index_type rstride[GFC_MAX_DIMENSIONS - 1];
+  index_type rstride0;
+  index_type roffset;
+  char *rptr;
+  char *dest;
+  /* s.* indicates the source array.  */
+  index_type sstride[GFC_MAX_DIMENSIONS - 1];
+  index_type sstride0;
+  index_type soffset;
+  const char *sptr;
+  const char *src;
+
+  index_type count[GFC_MAX_DIMENSIONS - 1];
+  index_type extent[GFC_MAX_DIMENSIONS - 1];
+  index_type dim;
+  index_type size;
+  index_type len;
+  index_type n;
+
+  if (!pbound)
+    pbound = zeros;
+
+  size = GFC_DESCRIPTOR_SIZE (ret);
+
+  which = which - 1;
+
+  extent[0] = 1;
+  count[0] = 0;
+  size = GFC_DESCRIPTOR_SIZE (array);
+  n = 0;
+  for (dim = 0; dim < GFC_DESCRIPTOR_RANK (array); dim++)
+    {
+      if (dim == which)
+        {
+          roffset = ret->dim[dim].stride * size;
+          if (roffset == 0)
+            roffset = size;
+          soffset = array->dim[dim].stride * size;
+          if (soffset == 0)
+            soffset = size;
+          len = array->dim[dim].ubound + 1 - array->dim[dim].lbound;
+        }
+      else
+        {
+          count[n] = 0;
+          extent[n] = array->dim[dim].ubound + 1 - array->dim[dim].lbound;
+          rstride[n] = ret->dim[dim].stride * size;
+          sstride[n] = array->dim[dim].stride * size;
+          n++;
+        }
+    }
+  if (sstride[0] == 0)
+    sstride[0] = size;
+  if (rstride[0] == 0)
+    rstride[0] = size;
+
+  dim = GFC_DESCRIPTOR_RANK (array);
+  rstride0 = rstride[0];
+  sstride0 = sstride[0];
+  rptr = ret->data;
+  sptr = array->data;
+  if (shift > 0)
+    len = len - shift;
+  else
+    len = len + shift;
+
+  while (rptr)
+    {
+      /* Do the shift for this dimension.  */
+      if (shift > 0)
+        {
+          src = &sptr[shift * soffset];
+          dest = rptr;
+        }
+      else
+        {
+          src = sptr;
+          dest = &rptr[-shift * roffset];
+        }
+      for (n = 0; n < len; n++)
+        {
+          memcpy (dest, src, size);
+          dest += roffset;
+          src += soffset;
+        }
+      if (shift >= 0)
+        {
+          n = shift;
+        }
+      else
+        {
+          dest = rptr;
+          n = -shift;
+        }
+
+      while (n--)
+        {
+          memcpy (dest, pbound, size);
+          dest += roffset;
+        }
+
+      /* Advance to the next section.  */
+      rptr += rstride0;
+      sptr += sstride0;
+      count[0]++;
+      n = 0;
+      while (count[n] == extent[n])
+        {
+          /* When we get to the end of a dimension, reset it and increment
+             the next dimension.  */
+          count[n] = 0;
+          /* We could precalculate these products, but this is a less
+             frequently used path so proabably not worth it.  */
+          rptr -= rstride[n] * extent[n];
+          sptr -= sstride[n] * extent[n];
+          n++;
+          if (n >= dim - 1)
+            {
+              /* Break out of the loop.  */
+              rptr = NULL;
+              break;
+            }
+          else
+            {
+              count[n]++;
+              rptr += rstride[n];
+              sptr += sstride[n];
+            }
+        }
+    }
+}
+
+
+void
+__eoshift0_4 (const gfc_array_char * ret, const gfc_array_char * array,
+    const GFC_INTEGER_4 * pshift, const char * pbound,
+    const GFC_INTEGER_4 * pdim)
+{
+  __eoshift0 (ret, array, *pshift, pbound, pdim ? *pdim : 1);
+}
+
+
+void
+__eoshift0_8 (const gfc_array_char * ret, const gfc_array_char * array,
+    const GFC_INTEGER_8 * pshift, const char * pbound,
+    const GFC_INTEGER_8 * pdim)
+{
+  __eoshift0 (ret, array, *pshift, pbound, pdim ? *pdim : 1);
+}
+
diff --git a/libgfortran/intrinsics/eoshift2.c b/libgfortran/intrinsics/eoshift2.c
new file mode 100644
index 00000000000..038588f78d2
--- /dev/null
+++ b/libgfortran/intrinsics/eoshift2.c
@@ -0,0 +1,204 @@
+/* Generic implementation of the RESHAPE intrinsic
+   Copyright 2002 Free Software Foundation, Inc.
+   Contributed by Paul Brook <paul@nowt.org>
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfor is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Ligbfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "config.h"
+#include <stdlib.h>
+#include <assert.h>
+#include <string.h>
+#include "libgfortran.h"
+
+static const char zeros[16] =
+  {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
+
+/* TODO: make this work for large shifts when
+   sizeof(int) < sizeof (index_type).  */
+
+static void
+__eoshift2 (const gfc_array_char * ret, const gfc_array_char * array,
+    int shift, const gfc_array_char * bound, int which)
+{
+  /* r.* indicates the return array.  */
+  index_type rstride[GFC_MAX_DIMENSIONS - 1];
+  index_type rstride0;
+  index_type roffset;
+  char *rptr;
+  char *dest;
+  /* s.* indicates the source array.  */
+  index_type sstride[GFC_MAX_DIMENSIONS - 1];
+  index_type sstride0;
+  index_type soffset;
+  const char *sptr;
+  const char *src;
+  /* b.* indicates the bound array.  */
+  index_type bstride[GFC_MAX_DIMENSIONS - 1];
+  index_type bstride0;
+  const char *bptr;
+
+  index_type count[GFC_MAX_DIMENSIONS - 1];
+  index_type extent[GFC_MAX_DIMENSIONS - 1];
+  index_type dim;
+  index_type size;
+  index_type len;
+  index_type n;
+
+  size = GFC_DESCRIPTOR_SIZE (ret);
+
+  which = which - 1;
+
+  extent[0] = 1;
+  count[0] = 0;
+  size = GFC_DESCRIPTOR_SIZE (array);
+  n = 0;
+  for (dim = 0; dim < GFC_DESCRIPTOR_RANK (array); dim++)
+    {
+      if (dim == which)
+        {
+          roffset = ret->dim[dim].stride * size;
+          if (roffset == 0)
+            roffset = size;
+          soffset = array->dim[dim].stride * size;
+          if (soffset == 0)
+            soffset = size;
+          len = array->dim[dim].ubound + 1 - array->dim[dim].lbound;
+        }
+      else
+        {
+          count[n] = 0;
+          extent[n] = array->dim[dim].ubound + 1 - array->dim[dim].lbound;
+          rstride[n] = ret->dim[dim].stride * size;
+          sstride[n] = array->dim[dim].stride * size;
+          if (bound)
+            bstride[n] = bound->dim[n].stride * size;
+          else
+            bstride[n] = 0;
+          n++;
+        }
+    }
+  if (sstride[0] == 0)
+    sstride[0] = size;
+  if (rstride[0] == 0)
+    rstride[0] = size;
+  if (bound && bstride[0] == 0)
+    bstride[0] = size;
+
+  dim = GFC_DESCRIPTOR_RANK (array);
+  rstride0 = rstride[0];
+  sstride0 = sstride[0];
+  bstride0 = bstride[0];
+  rptr = ret->data;
+  sptr = array->data;
+  if (bound)
+    bptr = bound->data;
+  else
+    bptr = zeros;
+
+  if (shift > 0)
+    len = len - shift;
+  else
+    len = len + shift;
+
+  while (rptr)
+    {
+      /* Do the shift for this dimension.  */
+      if (shift > 0)
+        {
+          src = &sptr[shift * soffset];
+          dest = rptr;
+        }
+      else
+        {
+          src = sptr;
+          dest = &rptr[-shift * roffset];
+        }
+      for (n = 0; n < len; n++)
+        {
+          memcpy (dest, src, size);
+          dest += roffset;
+          src += soffset;
+        }
+      if (shift >= 0)
+        {
+          n = shift;
+        }
+      else
+        {
+          dest = rptr;
+          n = -shift;
+        }
+
+      while (n--)
+        {
+          memcpy (dest, bptr, size);
+          dest += roffset;
+        }
+
+      /* Advance to the next section.  */
+      rptr += rstride0;
+      sptr += sstride0;
+      bptr += bstride0;
+      count[0]++;
+      n = 0;
+      while (count[n] == extent[n])
+        {
+          /* When we get to the end of a dimension, reset it and increment
+             the next dimension.  */
+          count[n] = 0;
+          /* We could precalculate these products, but this is a less
+             frequently used path so proabably not worth it.  */
+          rptr -= rstride[n] * extent[n];
+          sptr -= sstride[n] * extent[n];
+          bptr -= bstride[n] * extent[n];
+          n++;
+          if (n >= dim - 1)
+            {
+              /* Break out of the loop.  */
+              rptr = NULL;
+              break;
+            }
+          else
+            {
+              count[n]++;
+              rptr += rstride[n];
+              sptr += sstride[n];
+              bptr += bstride[n];
+            }
+        }
+    }
+}
+
+
+void
+__eoshift2_4 (const gfc_array_char * ret, const gfc_array_char * array,
+    const GFC_INTEGER_4 * pshift, const gfc_array_char * bound,
+    const GFC_INTEGER_4 * pdim)
+{
+  __eoshift2 (ret, array, *pshift, bound, pdim ? *pdim : 1);
+}
+
+
+void
+__eoshift2_8 (const gfc_array_char * ret, const gfc_array_char * array,
+    const GFC_INTEGER_8 * pshift, const gfc_array_char * bound,
+    const GFC_INTEGER_8 * pdim)
+{
+  __eoshift2 (ret, array, *pshift, bound, pdim ? *pdim : 1);
+}
+
diff --git a/libgfortran/intrinsics/ishftc.c b/libgfortran/intrinsics/ishftc.c
new file mode 100644
index 00000000000..0bb3d422eb2
--- /dev/null
+++ b/libgfortran/intrinsics/ishftc.c
@@ -0,0 +1,64 @@
+/* Implementation of ishftc intrinsic.
+   Copyright 2002 Free Software Foundation, Inc.
+   Contributed by Paul Brook <paul@nowt.org>
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfortran is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Libgfortran is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "libgfortran.h"
+
+#define ishftc4 prefix(ishftc4)
+GFC_INTEGER_4 ishftc4 (GFC_INTEGER_4, GFC_INTEGER_4, GFC_INTEGER_4);
+
+#define ishftc8 prefix(ishftc8)
+GFC_INTEGER_8 ishftc8 (GFC_INTEGER_8, GFC_INTEGER_8, GFC_INTEGER_8);
+
+GFC_INTEGER_4
+ishftc4 (GFC_INTEGER_4 i, GFC_INTEGER_4 shift, GFC_INTEGER_4 size)
+{
+  GFC_INTEGER_4 mask;
+  GFC_UINTEGER_4 bits;
+
+  if (shift < 0)
+    shift = shift + size;
+
+  if (shift == 0 || shift == size)
+    return i;
+
+  mask = (~(GFC_INTEGER_4)0) << size;
+  bits = i & ~mask;
+  return (i & mask) | (bits >> (size - shift)) | ((i << shift) & ~mask);
+}
+
+
+GFC_INTEGER_8
+ishftc8 (GFC_INTEGER_8 i, GFC_INTEGER_8 shift, GFC_INTEGER_8 size)
+{
+  GFC_INTEGER_8 mask;
+  GFC_UINTEGER_8 bits;
+
+  if (shift < 0)
+    shift = shift + size;
+
+  if (shift == 0 || shift == size)
+    return i;
+
+  mask = (~(GFC_INTEGER_8)0) << size;
+  bits = i & ~mask;
+  return (i & mask) | (bits >> (size - shift)) | ((i << shift) & ~mask);
+}
+
diff --git a/libgfortran/intrinsics/pack_generic.c b/libgfortran/intrinsics/pack_generic.c
new file mode 100644
index 00000000000..08c022e4e74
--- /dev/null
+++ b/libgfortran/intrinsics/pack_generic.c
@@ -0,0 +1,146 @@
+/* Generic implementation of the RESHAPE intrinsic
+   Copyright 2002 Free Software Foundation, Inc.
+   Contributed by Paul Brook <paul@nowt.org>
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfor is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Ligbfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "config.h"
+#include <stdlib.h>
+#include <assert.h>
+#include <string.h>
+#include "libgfortran.h"
+
+void
+__pack (const gfc_array_char * ret, const gfc_array_char * array,
+    const gfc_array_l4 * mask, const gfc_array_char * vector)
+{
+  /* r.* indicates the return array.  */
+  index_type rstride0;
+  char *rptr;
+  /* s.* indicates the source array.  */
+  index_type sstride[GFC_MAX_DIMENSIONS];
+  index_type sstride0;
+  const char *sptr;
+  /* m.* indicates the mask array.  */
+  index_type mstride[GFC_MAX_DIMENSIONS];
+  index_type mstride0;
+  const GFC_LOGICAL_4 *mptr;
+
+  index_type count[GFC_MAX_DIMENSIONS];
+  index_type extent[GFC_MAX_DIMENSIONS];
+  index_type n;
+  index_type dim;
+  index_type size;
+  index_type nelem;
+
+  size = GFC_DESCRIPTOR_SIZE (array);
+  dim = GFC_DESCRIPTOR_RANK (array);
+  for (n = 0; n < dim; n++)
+    {
+      count[n] = 0;
+      extent[n] = array->dim[n].ubound + 1 - array->dim[n].lbound;
+      sstride[n] = array->dim[n].stride * size;
+      mstride[n] = mask->dim[n].stride;
+    }
+  if (sstride[0] == 0)
+    sstride[0] = size;
+  if (mstride[0] == 0)
+    mstride[0] = 1;
+
+  rstride0 = ret->dim[0].stride * size;
+  if (rstride0 == 0)
+    rstride0 = size;
+  sstride0 = sstride[0];
+  mstride0 = mstride[0];
+  rptr = ret->data;
+  sptr = array->data;
+  mptr = mask->data;
+
+  /* Use the same loop for both logical types. */
+  if (GFC_DESCRIPTOR_SIZE (mask) != 4)
+    {
+      if (GFC_DESCRIPTOR_SIZE (mask) != 8)
+        runtime_error ("Funny sized logical array");
+      for (n = 0; n < dim; n++)
+        mstride[n] <<= 1;
+      mstride0 <<= 1;
+      mptr = GFOR_POINTER_L8_TO_L4 (mptr);
+    }
+
+  while (sptr)
+    {
+      /* Test this element.  */
+      if (*mptr)
+        {
+          /* Add it.  */
+          memcpy (rptr, sptr, size);
+          rptr += rstride0;
+        }
+      /* Advance to the next element.  */
+      sptr += sstride0;
+      mptr += mstride0;
+      count[0]++;
+      n = 0;
+      while (count[n] == extent[n])
+        {
+          /* When we get to the end of a dimension, reset it and increment
+             the next dimension.  */
+          count[n] = 0;
+          /* We could precalculate these products, but this is a less
+             frequently used path so proabably not worth it.  */
+          sptr -= sstride[n] * extent[n];
+          mptr -= mstride[n] * extent[n];
+          n++;
+          if (n >= dim)
+            {
+              /* Break out of the loop.  */
+              sptr = NULL;
+              break;
+            }
+          else
+            {
+              count[n]++;
+              sptr += sstride[n];
+              mptr += mstride[n];
+            }
+        }
+    }
+
+  /* Add any remaining elements from VECTOR.  */
+  if (vector)
+    {
+      n = vector->dim[0].ubound + 1 - vector->dim[0].lbound;
+      nelem = ((rptr - ret->data) / rstride0);
+      if (n > nelem)
+        {
+          sstride0 = vector->dim[0].stride * size;
+          if (sstride0 == 0)
+            sstride0 = size;
+
+          sptr = vector->data + sstride0 * nelem;
+          n -= nelem;
+          while (n--)
+            {
+              memcpy (rptr, sptr, size);
+              rptr += rstride0;
+              sptr += sstride0;
+            }
+        }
+    }
+}
+
diff --git a/libgfortran/intrinsics/random.c b/libgfortran/intrinsics/random.c
new file mode 100644
index 00000000000..b578148f469
--- /dev/null
+++ b/libgfortran/intrinsics/random.c
@@ -0,0 +1,362 @@
+/* Implementation of the RANDOM intrinsics
+   Copyright 2002 Free Software Foundation, Inc.
+   Contributed by Lars Segerlund <seger@linuxmail.org>
+
+  The algorithm was taken from the paper :
+
+	Mersenne Twister:	623-dimensionally equidistributed
+				uniform pseudorandom generator.
+
+	by:	Makoto Matsumoto
+		Takuji Nishimura
+
+	Which appeared in the: ACM Transactions on Modelling and Computer
+	Simulations: Special Issue on Uniform Random Number
+	Generation. ( Early in 1998 ).
+
+This file is part of the GNU Fortran 95 runtime library (libgfortran).
+
+Libgfortran is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Ligbfortran is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "config.h"
+#include <stdio.h>
+#include <stdlib.h>
+#include <sys/types.h>
+#include <sys/stat.h>
+#include <fcntl.h>
+#include <assert.h>
+#include "libgfortran.h"
+
+/*Use the 'big' generator by default ( period -> 2**19937 ).  */
+
+#define MT19937
+
+/* Define the necessary constants for the algorithm.  */
+
+#ifdef  MT19937
+enum constants
+{
+  N = 624, M = 397, R = 19, TU = 11, TS = 7, TT = 15, TL = 17
+};
+#define M_A	0x9908B0DF
+#define T_B	0x9D2C5680
+#define T_C	0xEFC60000
+#else
+enum constants
+{
+  N = 351, M = 175, R = 19, TU = 11, TS = 7, TT = 15, TL = 17
+};
+#define M_A	0xE4BD75F5
+#define T_B	0x655E5280
+#define T_C	0xFFD58000
+#endif
+
+static int i = N;
+static unsigned int seed[N];
+
+/* This is the routine which handles the seeding of the generator,
+   and also reading and writing of the seed.  */
+
+#define random_seed prefix(random_seed)
+void
+random_seed (GFC_INTEGER_4 * size, const gfc_array_i4 * put,
+	     const gfc_array_i4 * get)
+{
+  /* Initialize the seed in system dependent manner.  */
+  if (get == NULL && put == NULL && size == NULL)
+    {
+      int fd;
+      fd = open ("/dev/urandom", O_RDONLY);
+      if (fd == 0)
+	{
+	  /* We dont have urandom.  */
+	  GFC_UINTEGER_4 s = (GFC_UINTEGER_4) seed;
+	  for (i = 0; i < N; i++)
+	    {
+	      s = s * 29943829 - 1;
+	      seed[i] = s;
+	    }
+	}
+      else
+	{
+	  /* Using urandom, might have a length issue.  */
+	  read (fd, &seed[0], sizeof (GFC_UINTEGER_4) * N);
+	  close (fd);
+	}
+      return;
+    }
+
+  /* Return the size of the seed */
+  if (size != NULL)
+    {
+      *size = N;
+      return;
+    }
+
+  /* if we have gotten to this pount we have a get or put
+   * now we check it the array fulfills the demands in the standard .
+   */
+
+  /* Set the seed to PUT data */
+  if (put != NULL)
+    {
+      /* if the rank of the array is not 1 abort */
+      if (GFC_DESCRIPTOR_RANK (put) != 1)
+	abort ();
+
+      /* if the array is too small abort */
+      if (((put->dim[0].ubound + 1 - put->dim[0].lbound)) < N)
+	abort ();
+
+      /* If this is the case the array is a temporary */
+      if (get->dim[0].stride == 0)
+	return;
+
+      /*  This code now should do correct strides. */
+      for (i = 0; i < N; i++)
+	seed[i] = put->data[i * put->dim[0].stride];
+    }
+
+  /* Return the seed to GET data */
+  if (get != NULL)
+    {
+      /* if the rank of the array is not 1 abort */
+      if (GFC_DESCRIPTOR_RANK (get) != 1)
+	abort ();
+
+      /* if the array is too small abort */
+      if (((get->dim[0].ubound + 1 - get->dim[0].lbound)) < N)
+	abort ();
+
+      /* If this is the case the array is a temporary */
+      if (get->dim[0].stride == 0)
+	return;
+
+      /*  This code now should do correct strides. */
+      for (i = 0; i < N; i++)
+	get->data[i * get->dim[0].stride] = seed[i];
+    }
+}
+
+/* Here is the internal routine which generates the random numbers
+   in 'batches' based upon the need for a new batch.
+   It's an integer based routine known as 'Mersenne Twister'.
+   This implementation still lacks 'tempering' and a good verification,
+   but gives very good metrics.  */
+
+static void
+random_generate (void)
+{
+  /* 32 bits.  */
+  GFC_UINTEGER_4 y;
+
+  /* Generate batch of N.  */
+  int k, m;
+  for (k = 0, m = M; k < N - 1; k++)
+    {
+      y = (seed[k] & (-1 << R)) | (seed[k + 1] & ((1u << R) - 1));
+      seed[k] = seed[m] ^ (y >> 1) ^ (-(GFC_INTEGER_4) (y & 1) & M_A);
+      if (++m >= N)
+	m = 0;
+    }
+
+  y = (seed[N - 1] & (-1 << R)) | (seed[0] & ((1u << R) - 1));
+  seed[N - 1] = seed[M - 1] ^ (y >> 1) ^ (-(GFC_INTEGER_4) (y & 1) & M_A);
+  i = 0;
+}
+
+/* A routine to return a REAL(KIND=4).  */
+
+#define random_r4 prefix(random_r4)
+void
+random_r4 (GFC_REAL_4 * harv)
+{
+  /* Regenerate if we need to.  */
+  if (i >= N)
+    random_generate ();
+
+  /* Convert uint32 to REAL(KIND=4).  */
+  *harv = (GFC_REAL_4) ((GFC_REAL_4) (GFC_UINTEGER_4) seed[i++] /
+			(GFC_REAL_4) (~(GFC_UINTEGER_4) 0));
+}
+
+/* A routine to return a REAL(KIND=8).  */
+
+#define random_r8 prefix(random_r8)
+void
+random_r8 (GFC_REAL_8 * harv)
+{
+  /* Regenerate if we need to, may waste one 32-bit value.  */
+  if ((i + 1) >= N)
+    random_generate ();
+
+  /* Convert two uint32 to a REAL(KIND=8).  */
+  *harv = ((GFC_REAL_8) ((((GFC_UINTEGER_8) seed[i+1]) << 32) + seed[i])) /
+	  (GFC_REAL_8) (~(GFC_UINTEGER_8) 0);
+  i += 2;
+}
+
+/* Code to handle arrays will follow here.  */
+
+/* REAL(KIND=4) REAL array.  */
+
+#define arandom_r4 prefix(arandom_r4)
+void
+arandom_r4 (gfc_array_r4 * harv)
+{
+  index_type count[GFC_MAX_DIMENSIONS - 1];
+  index_type extent[GFC_MAX_DIMENSIONS - 1];
+  index_type stride[GFC_MAX_DIMENSIONS - 1];
+  index_type stride0;
+  index_type dim;
+  GFC_REAL_4 *dest;
+  int n;
+
+  dest = harv->data;
+
+  if (harv->dim[0].stride == 0)
+    harv->dim[0].stride = 1;
+
+  dim = GFC_DESCRIPTOR_RANK (harv);
+
+  for (n = 0; n < dim; n++)
+    {
+      count[n] = 0;
+      stride[n] = harv->dim[n].stride;
+      extent[n] = harv->dim[n].ubound + 1 - harv->dim[n].lbound;
+      if (extent[n] <= 0)
+	return;
+    }
+
+  stride0 = stride[0];
+
+  while (dest)
+    {
+      /* Set the elements.  */
+
+      /* regenerate if we need to */
+      if (i >= N)
+	random_generate ();
+
+      /* Convert uint32 to float in a hopefully g95 compiant manner */
+      *dest = (GFC_REAL_4) ((GFC_REAL_4) (GFC_UINTEGER_4) seed[i++] /
+			    (GFC_REAL_4) (~(GFC_UINTEGER_4) 0));
+
+      /* Advance to the next element.  */
+      dest += stride0;
+      count[0]++;
+      /* Advance to the next source element.  */
+      n = 0;
+      while (count[n] == extent[n])
+	{
+	  /* When we get to the end of a dimension,
+	     reset it and increment
+	     the next dimension.  */
+	  count[n] = 0;
+	  /* We could precalculate these products,
+	     but this is a less
+	     frequently used path so proabably not worth it.  */
+	  dest -= stride[n] * extent[n];
+	  n++;
+	  if (n == dim)
+	    {
+	      dest = NULL;
+	      break;
+	    }
+	  else
+	    {
+	      count[n]++;
+	      dest += stride[n];
+	    }
+	}
+    }
+}
+
+/* REAL(KIND=8) array.  */
+
+#define arandom_r8 prefix(arandom_r8)
+void
+arandom_r8 (gfc_array_r8 * harv)
+{
+  index_type count[GFC_MAX_DIMENSIONS - 1];
+  index_type extent[GFC_MAX_DIMENSIONS - 1];
+  index_type stride[GFC_MAX_DIMENSIONS - 1];
+  index_type stride0;
+  index_type dim;
+  GFC_REAL_8 *dest;
+  int n;
+
+  dest = harv->data;
+
+  if (harv->dim[0].stride == 0)
+    harv->dim[0].stride = 1;
+
+  dim = GFC_DESCRIPTOR_RANK (harv);
+
+  for (n = 0; n < dim; n++)
+    {
+      count[n] = 0;
+      stride[n] = harv->dim[n].stride;
+      extent[n] = harv->dim[n].ubound + 1 - harv->dim[n].lbound;
+      if (extent[n] <= 0)
+	return;
+    }
+
+  stride0 = stride[0];
+
+  while (dest)
+    {
+      /* Set the elements.  */
+
+      /* regenerate if we need to, may waste one 32-bit value */
+      if ((i + 1) >= N)
+	random_generate ();
+
+      /* Convert two uint32 to a REAL(KIND=8).  */
+      *dest = ((GFC_REAL_8) ((((GFC_UINTEGER_8) seed[i+1]) << 32) + seed[i])) /
+	      (GFC_REAL_8) (~(GFC_UINTEGER_8) 0);
+      i += 2;
+
+      /* Advance to the next element.  */
+      dest += stride0;
+      count[0]++;
+      /* Advance to the next source element.  */
+      n = 0;
+      while (count[n] == extent[n])
+	{
+	  /* When we get to the end of a dimension,
+	     reset it and increment
+	     the next dimension.  */
+	  count[n] = 0;
+	  /* We could precalculate these products,
+	     but this is a less
+	     frequently used path so proabably not worth it.  */
+	  dest -= stride[n] * extent[n];
+	  n++;
+	  if (n == dim)
+	    {
+	      dest = NULL;
+	      break;
+	    }
+	  else
+	    {
+	      count[n]++;
+	      dest += stride[n];
+	    }
+	}
+    }
+}
+
diff --git a/libgfortran/intrinsics/reshape_generic.c b/libgfortran/intrinsics/reshape_generic.c
new file mode 100644
index 00000000000..ca6f6aacd00
--- /dev/null
+++ b/libgfortran/intrinsics/reshape_generic.c
@@ -0,0 +1,231 @@
+/* Generic implementation of the RESHAPE intrinsic
+   Copyright 2002 Free Software Foundation, Inc.
+   Contributed by Paul Brook <paul@nowt.org>
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfor is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Ligbfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "config.h"
+#include <stdlib.h>
+#include <string.h>
+#include <assert.h>
+#include "libgfortran.h"
+
+typedef GFC_ARRAY_DESCRIPTOR(1, index_type) shape_type;
+typedef GFC_ARRAY_DESCRIPTOR(GFC_MAX_DIMENSIONS, char) parray;
+
+
+/* The shape parameter is ignored. We can currently deduce the shape from the
+   return array.  */
+
+void
+__reshape (parray * ret, parray * source, shape_type * shape,
+           parray * pad, shape_type * order)
+{
+  /* r.* indicates the return array.  */
+  index_type rcount[GFC_MAX_DIMENSIONS - 1];
+  index_type rextent[GFC_MAX_DIMENSIONS - 1];
+  index_type rstride[GFC_MAX_DIMENSIONS - 1];
+  index_type rstride0;
+  index_type rdim;
+  index_type rsize;
+  char *rptr;
+  /* s.* indicates the source array.  */
+  index_type scount[GFC_MAX_DIMENSIONS - 1];
+  index_type sextent[GFC_MAX_DIMENSIONS - 1];
+  index_type sstride[GFC_MAX_DIMENSIONS - 1];
+  index_type sstride0;
+  index_type sdim;
+  index_type ssize;
+  const char *sptr;
+  /* p.* indicates the pad array.  */
+  index_type pcount[GFC_MAX_DIMENSIONS - 1];
+  index_type pextent[GFC_MAX_DIMENSIONS - 1];
+  index_type pstride[GFC_MAX_DIMENSIONS - 1];
+  index_type pdim;
+  index_type psize;
+  const char *pptr;
+
+  const char *src;
+  int n;
+  int dim;
+  int size;
+
+  size = GFC_DESCRIPTOR_SIZE (ret);
+  if (ret->dim[0].stride == 0)
+    ret->dim[0].stride = 1;
+  if (source->dim[0].stride == 0)
+    source->dim[0].stride = 1;
+  if (shape->dim[0].stride == 0)
+    shape->dim[0].stride = 1;
+  if (pad && pad->dim[0].stride == 0)
+    pad->dim[0].stride = 1;
+  if (order && order->dim[0].stride == 0)
+    order->dim[0].stride = 1;
+
+  rdim = GFC_DESCRIPTOR_RANK (ret);
+  rsize = 1;
+  for (n = 0; n < rdim; n++)
+    {
+      if (order)
+        dim = order->data[n * order->dim[0].stride] - 1;
+      else
+        dim = n;
+
+      rcount[n] = 0;
+      rstride[n] = ret->dim[dim].stride;
+      rextent[n] = ret->dim[dim].ubound + 1 - ret->dim[dim].lbound;
+
+      if (rextent[n] != shape->data[dim * shape->dim[0].stride])
+        runtime_error ("shape and target do not conform");
+
+      if (rsize == rstride[n])
+        rsize *= rextent[n];
+      else
+        rsize = 0;
+      if (rextent[dim] <= 0)
+        return;
+    }
+
+  sdim = GFC_DESCRIPTOR_RANK (source);
+  ssize = 1;
+  for (n = 0; n < sdim; n++)
+    {
+      scount[n] = 0;
+      sstride[n] = source->dim[n].stride;
+      sextent[n] = source->dim[n].ubound + 1 - source->dim[n].lbound;
+      if (sextent[n] <= 0)
+        abort ();
+
+      if (rsize == sstride[n])
+        ssize *= sextent[n];
+      else
+        ssize = 0;
+    }
+
+  if (pad)
+    {
+      if (pad->dim[0].stride == 0)
+        pad->dim[0].stride = 1;
+      pdim = GFC_DESCRIPTOR_RANK (pad);
+      psize = 1;
+      for (n = 0; n < pdim; n++)
+        {
+          pcount[n] = 0;
+          pstride[n] = pad->dim[n].stride;
+          pextent[n] = pad->dim[n].ubound + 1 - pad->dim[n].lbound;
+          if (pextent[n] <= 0)
+            abort ();
+          if (psize == pstride[n])
+            psize *= pextent[n];
+          else
+            rsize = 0;
+        }
+      pptr = pad->data;
+    }
+  else
+    {
+      pdim = 0;
+      psize = 1;
+      pptr = NULL;
+    }
+
+  if (rsize != 0 && ssize != 0 && psize != 0)
+    {
+      rsize *= size;
+      ssize *= size;
+      psize *= size;
+      reshape_packed (ret->data, rsize, source->data, ssize,
+		      pad ? pad->data : NULL, psize);
+      return;
+    }
+  rptr = ret->data;
+  src = sptr = source->data;
+  rstride0 = rstride[0] * size;
+  sstride0 = sstride[0] * size;
+
+  while (rptr)
+    {
+      /* Select between the source and pad arrays.  */
+      memcpy(rptr, src, size);
+      /* Advance to the next element.  */
+      rptr += rstride0;
+      src += sstride0;
+      rcount[0]++;
+      scount[0]++;
+      /* Advance to the next destination element.  */
+      n = 0;
+      while (rcount[n] == rextent[n])
+        {
+          /* When we get to the end of a dimension, reset it and increment
+             the next dimension.  */
+          rcount[n] = 0;
+          /* We could precalculate these products, but this is a less
+             frequently used path so proabably not worth it.  */
+          rptr -= rstride[n] * rextent[n] * size;
+          n++;
+          if (n == rdim)
+            {
+              /* Break out of the loop.  */
+              rptr = NULL;
+              break;
+            }
+          else
+            {
+              rcount[n]++;
+              rptr += rstride[n] * size;
+            }
+        }
+      /* Advance to the next source element.  */
+      n = 0;
+      while (scount[n] == sextent[n])
+        {
+          /* When we get to the end of a dimension, reset it and increment
+             the next dimension.  */
+          scount[n] = 0;
+          /* We could precalculate these products, but this is a less
+             frequently used path so proabably not worth it.  */
+          src -= sstride[n] * sextent[n] * size;
+          n++;
+          if (n == sdim)
+            {
+              if (sptr && pad)
+                {
+                  /* Switch to the pad array.  */
+                  sptr = NULL;
+                  sdim = pdim;
+                  for (dim = 0; dim < pdim; dim++)
+                    {
+                      scount[dim] = pcount[dim];
+                      sextent[dim] = pextent[dim];
+                      sstride[dim] = pstride[dim];
+                      sstride0 = sstride[0] * size;
+                    }
+                }
+              /* We now start again from the beginning of the pad array.  */
+              src = pptr;
+              break;
+            }
+          else
+            {
+              scount[n]++;
+              sptr += sstride[n] * size;
+            }
+        }
+    }
+}
+
diff --git a/libgfortran/intrinsics/reshape_packed.c b/libgfortran/intrinsics/reshape_packed.c
new file mode 100644
index 00000000000..eef885feb85
--- /dev/null
+++ b/libgfortran/intrinsics/reshape_packed.c
@@ -0,0 +1,46 @@
+/* Implementation of the RESHAPE intrinsic for packed arrays
+   Copyright 2002 Free Software Foundation, Inc.
+   Contributed by Paul Brook <paul@nowt.org>
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfor is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Ligbfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "config.h"
+#include "libgfortran.h"
+
+#include <string.h>
+
+/* Reshape function where all arrays are packed.  Basically just memcpy.  */
+
+void
+reshape_packed (char * ret, index_type rsize, const char * source,
+		index_type ssize, const char * pad, index_type psize)
+{
+  index_type size;
+
+  size = (rsize > ssize) ? ssize : rsize;
+  memcpy (ret, source, size);
+  ret += size;
+  rsize -= size;
+  while (rsize > 0)
+    {
+      size = (rsize > psize) ? psize : rsize;
+      memcpy (ret, pad, size);
+      ret += size;
+      rsize -= size;
+    }
+}
diff --git a/libgfortran/intrinsics/selected_kind.f90 b/libgfortran/intrinsics/selected_kind.f90
new file mode 100644
index 00000000000..62d11c7f596
--- /dev/null
+++ b/libgfortran/intrinsics/selected_kind.f90
@@ -0,0 +1,90 @@
+!   Copyright 2003 Free Software Foundation, Inc.
+!   Contributed by Kejia Zhao <kejia_zh@yahoo.com.cn>
+!
+!This file is part of the GNU Fortran 95 runtime library (libgfor).
+!
+!GNU libgfor is free software; you can redistribute it and/or
+!modify it under the terms of the GNU Lesser General Public
+!License as published by the Free Software Foundation; either
+!version 2.1 of the License, or (at your option) any later version.
+!
+!GNU libgfor is distributed in the hope that it will be useful,
+!but WITHOUT ANY WARRANTY; without even the implied warranty of
+!MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+!GNU Lesser General Public License for more details.
+!
+!You should have received a copy of the GNU Lesser General Public
+!License along with libgfor; see the file COPYING.  If not,
+!write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+!Boston, MA 02111-1307, USA.
+!
+
+function selected_int_kind (r)
+  implicit none
+  integer, intent (in) :: r
+  integer :: selected_int_kind
+  integer :: i
+  ! Integer kind_range table
+  integer, parameter :: c = 4
+  type :: int_info
+    integer :: kind
+    integer :: range
+  end type int_info
+  type (int_info), parameter :: int_infos (c) = &
+      (/int_info (1, range (0_1)), &
+        int_info (2, range (0_2)), &
+        int_info (4, range (0_4)), &
+        int_info (8, range (0_8))/)
+
+  do i = 1, c
+    if (r <= int_infos (i) % range) then
+      selected_int_kind = int_infos (i) % kind
+      return
+    end if
+  end do
+  selected_int_kind = -1
+  return
+end function
+
+function selected_real_kind (p, r)
+  implicit none
+  integer, optional, intent (in) :: p, r
+  integer :: selected_real_kind
+  integer :: i, p2, r2
+  logical :: found_p, found_r
+  ! Real kind_precision_range table
+  integer, parameter :: c = 2
+  type :: real_info
+    integer :: kind
+    integer :: precision
+    integer :: range
+  end type real_info
+  type (real_info) :: real_infos (c) = &
+      (/real_info (4, precision (0.0_4), range (0.0_4)), &
+        real_info (8, precision (0.0_8), range (0.0_8))/)
+
+  selected_real_kind = 0
+  p2 = 0
+  r2 = 0
+  found_p = .false.
+  found_r = .false.
+
+  if (present (p)) p2 = p
+  if (present (r)) r2 = r
+
+  ! Assumes each type has a greater precision and range than previous one.
+
+  do i = 1, c
+    if (p2 <= real_infos (i) % precision) found_p = .true.
+    if (r2 <= real_infos (i) % range) found_r = .true.
+    if (found_p .and. found_r) then
+      selected_real_kind = real_infos (i) % kind
+      return
+    end if
+  end do
+
+  if (.not. (found_p)) selected_real_kind = selected_real_kind - 1
+  if (.not. (found_r)) selected_real_kind = selected_real_kind - 2
+
+  return
+end function
diff --git a/libgfortran/intrinsics/size.c b/libgfortran/intrinsics/size.c
new file mode 100644
index 00000000000..5664c304988
--- /dev/null
+++ b/libgfortran/intrinsics/size.c
@@ -0,0 +1,56 @@
+/* Implementation of the size intrinsic.
+   Copyright 2002 Free Software Foundation, Inc.
+   Contributed by Paul Brook <paul@nowt.org>
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfor is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Libgfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "libgfortran.h"
+
+typedef GFC_ARRAY_DESCRIPTOR (GFC_MAX_DIMENSIONS, void) array_t;
+
+#define size0 prefix(size0)
+index_type size0 (const array_t * array)
+{
+  int n;
+  index_type size;
+  index_type len;
+
+  size = 1;
+  for (n = 0; n < GFC_DESCRIPTOR_RANK (array); n++)
+    {
+      len = array->dim[n].ubound + 1 - array->dim[n].lbound;
+      if (len < 0)
+        len = 0;
+      size *= len;
+    }
+  return size;
+}
+
+#define size1 prefix(size1)
+index_type size1 (const array_t * array, index_type dim)
+{
+  index_type size;
+
+  dim--;
+
+  size = array->dim[dim].ubound + 1 - array->dim[dim].lbound;
+  if (size < 0)
+    size = 0;
+  return size;
+}
+
diff --git a/libgfortran/intrinsics/spread_generic.c b/libgfortran/intrinsics/spread_generic.c
new file mode 100644
index 00000000000..a789c98448b
--- /dev/null
+++ b/libgfortran/intrinsics/spread_generic.c
@@ -0,0 +1,118 @@
+/* Generic implementation of the RESHAPE intrinsic
+   Copyright 2002 Free Software Foundation, Inc.
+   Contributed by Paul Brook <paul@nowt.org>
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfor is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Ligbfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "config.h"
+#include <stdlib.h>
+#include <assert.h>
+#include <string.h>
+#include "libgfortran.h"
+
+void
+__spread (const gfc_array_char * ret, const gfc_array_char * source,
+    const index_type * along, const index_type * pncopies)
+{
+  /* r.* indicates the return array.  */
+  index_type rstride[GFC_MAX_DIMENSIONS - 1];
+  index_type rstride0;
+  index_type rdelta;
+  char *rptr;
+  char *dest;
+  /* s.* indicates the source array.  */
+  index_type sstride[GFC_MAX_DIMENSIONS - 1];
+  index_type sstride0;
+  const char *sptr;
+
+  index_type count[GFC_MAX_DIMENSIONS - 1];
+  index_type extent[GFC_MAX_DIMENSIONS - 1];
+  index_type n;
+  index_type dim;
+  index_type size;
+  index_type ncopies;
+
+  size = GFC_DESCRIPTOR_SIZE (source);
+  dim = 0;
+  for (n = 0; n < GFC_DESCRIPTOR_RANK (ret); n++)
+    {
+      if (n == *along - 1)
+        {
+          rdelta = ret->dim[n].stride * size;
+        }
+      else
+        {
+          count[dim] = 0;
+          extent[dim] = source->dim[dim].ubound + 1 - source->dim[dim].lbound;
+          sstride[dim] = source->dim[dim].stride * size;
+          rstride[dim] = ret->dim[n].stride * size;
+          dim++;
+        }
+    }
+  dim = GFC_DESCRIPTOR_RANK (source);
+  if (sstride[0] == 0)
+    sstride[0] = size;
+  if (rstride[0] == 0)
+    rstride[0] = size;
+
+  sstride0 = sstride[0];
+  rstride0 = rstride[0];
+  rptr = ret->data;
+  sptr = source->data;
+  ncopies = *pncopies;
+
+  while (sptr)
+    {
+      /* Spread this element.  */
+      dest = rptr;
+      for (n = 0; n < ncopies; n++)
+        {
+          memcpy (dest, sptr, size);
+          dest += rdelta;
+        }
+      /* Advance to the next element.  */
+      sptr += sstride0;
+      rptr += rstride0;
+      count[0]++;
+      n = 0;
+      while (count[n] == extent[n])
+        {
+          /* When we get to the end of a dimension, reset it and increment
+             the next dimension.  */
+          count[n] = 0;
+          /* We could precalculate these products, but this is a less
+             frequently used path so probably not worth it.  */
+          sptr -= sstride[n] * extent[n];
+          rptr -= rstride[n] * extent[n];
+          n++;
+          if (n >= dim)
+            {
+              /* Break out of the loop.  */
+              sptr = NULL;
+              break;
+            }
+          else
+            {
+              count[n]++;
+              sptr += sstride[n];
+              rptr += rstride[n];
+            }
+        }
+    }
+}
+
diff --git a/libgfortran/intrinsics/string_intrinsics.c b/libgfortran/intrinsics/string_intrinsics.c
new file mode 100644
index 00000000000..999807ed1d0
--- /dev/null
+++ b/libgfortran/intrinsics/string_intrinsics.c
@@ -0,0 +1,394 @@
+/* String intrinsics helper functions.
+   Copyright 2002 Free Software Foundation, Inc.
+   Contributed by Paul Brook <paul@nowt.org>
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfortran is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Libgfortran is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+
+/* Unlike what the name of this file suggests, we don't actually
+   implement the Fortran intrinsics here.  At least, not with the
+   names they have in the standard.  The functions here provide all
+   the support we need for the standard string intrinsics, and the
+   compiler translates the actual intrinsics calls to calls to
+   functions in this file.  */
+
+#include <stdlib.h>
+#include <string.h>
+
+#include "libgfortran.h"
+
+
+/* String functions.  */
+
+#define copy_string prefix(copy_string)
+void copy_string (GFC_INTEGER_4, char *, GFC_INTEGER_4, const char *);
+
+#define concat_string prefix(concat_string)
+void concat_string (GFC_INTEGER_4, char *,
+		    GFC_INTEGER_4, const char *,
+		    GFC_INTEGER_4, const char *);
+
+#define string_len_trim prefix(string_len_trim)
+GFC_INTEGER_4 string_len_trim (GFC_INTEGER_4, const char *);
+
+#define adjustl prefix(adjustl)
+void adjustl (char *, GFC_INTEGER_4, const char *);
+
+#define adjustr prefix(adjustr)
+void adjustr (char *, GFC_INTEGER_4, const char *);
+
+#define string_index prefix(string_index)
+GFC_INTEGER_4 string_index (GFC_INTEGER_4, const char *, GFC_INTEGER_4,
+			    const char *, GFC_LOGICAL_4);
+
+#define string_scan prefix(string_scan)
+GFC_INTEGER_4 string_scan (GFC_INTEGER_4, const char *, GFC_INTEGER_4,
+                           const char *, GFC_LOGICAL_4);
+
+#define string_verify prefix(string_verify)
+GFC_INTEGER_4 string_verify (GFC_INTEGER_4, const char *, GFC_INTEGER_4,
+                             const char *, GFC_LOGICAL_4);
+
+#define string_trim prefix(string_trim)
+void string_trim (GFC_INTEGER_4 *, void **, GFC_INTEGER_4, const char *);
+
+#define string_repeat prefix(string_repeat)
+void string_repeat (char *, GFC_INTEGER_4, const char *, GFC_INTEGER_4);
+
+/* The two areas may overlap so we use memmove.  */
+
+void
+copy_string (GFC_INTEGER_4 destlen, char * dest,
+	     GFC_INTEGER_4 srclen, const char * src)
+{
+  if (srclen >= destlen)
+    {
+      /* This will truncate if too long.  */
+      memmove (dest, src, destlen);
+      /*memcpy (dest, src, destlen);*/
+    }
+  else
+    {
+      memmove (dest, src, srclen);
+      /*memcpy (dest, src, srclen);*/
+      /* Pad with spaces.  */
+      memset (&dest[srclen], ' ', destlen - srclen);
+    }
+}
+
+
+/* Strings of unequal length are extended with pad characters.  */
+
+GFC_INTEGER_4
+compare_string (GFC_INTEGER_4 len1, const char * s1,
+		GFC_INTEGER_4 len2, const char * s2)
+{
+  int res;
+  const char *s;
+  int len;
+
+  res = strncmp (s1, s2, (len1 < len2) ? len1 : len2);
+  if (res != 0)
+    return res;
+
+  if (len1 == len2)
+    return 0;
+
+  if (len1 < len2)
+    {
+      len = len2 - len1;
+      s = &s2[len1];
+      res = -1;
+    }
+  else
+    {
+      len = len1 - len2;
+      s = &s1[len2];
+      res = 1;
+    }
+
+  while (len--)
+    {
+      if (*s != ' ')
+        {
+          if (*s > ' ')
+            return res;
+          else
+            return -res;
+        }
+      s++;
+    }
+
+  return 0;
+}
+
+
+/* The destination and source should not overlap.  */
+
+void
+concat_string (GFC_INTEGER_4 destlen, char * dest,
+	       GFC_INTEGER_4 len1, const char * s1,
+	       GFC_INTEGER_4 len2, const char * s2)
+{
+  if (len1 >= destlen)
+    {
+      memcpy (dest, s1, destlen);
+      return;
+    }
+  memcpy (dest, s1, len1);
+  dest += len1;
+  destlen -= len1;
+
+  if (len2 >= destlen)
+    {
+      memcpy (dest, s2, destlen);
+      return;
+    }
+
+  memcpy (dest, s2, len2);
+  memset (&dest[len2], ' ', destlen - len2);
+}
+
+
+/* Return string with all trailing blanks removed.  */
+
+void
+string_trim (GFC_INTEGER_4 * len, void ** dest, GFC_INTEGER_4 slen, const char * src)
+{
+  int i;
+
+  /* Determine length of result string.  */
+  for (i = slen - 1; i >= 0; i--)
+    {
+      if (src[i] != ' ')
+        break;
+    }
+  *len = i + 1;
+
+  if (*len > 0)
+    {
+      /* Allocate space for result string.  */
+      *dest = internal_malloc (*len);
+
+      /* copy string if necessary.  */
+      memmove (*dest, src, *len);
+    }
+}
+
+
+/* The length of a string not including trailing blanks.  */
+
+GFC_INTEGER_4
+string_len_trim (GFC_INTEGER_4 len, const char * s)
+{
+  int i;
+
+  for (i = len - 1; i >= 0; i--)
+    {
+      if (s[i] != ' ')
+        break;
+    }
+  return i + 1;
+}
+
+
+/* Find a substring within a string.  */
+
+GFC_INTEGER_4
+string_index (GFC_INTEGER_4 slen, const char * str, GFC_INTEGER_4 sslen,
+	      const char * sstr, GFC_LOGICAL_4 back)
+{
+  int start;
+  int last;
+  int i;
+  int delta;
+
+  if (sslen == 0)
+    return 1;
+
+  if (!back)
+    {
+      last = slen + 1 - sslen;
+      start = 0;
+      delta = 1;
+    }
+  else
+    {
+      last = -1;
+      start = slen - sslen;
+      delta = -1;
+    }
+  i = 0;
+  for (; start != last; start+= delta)
+    {
+      for (i = 0; i < sslen; i++)
+        {
+          if (str[start + i] != sstr[i])
+            break;
+        }
+      if (i == sslen)
+        return (start + 1);
+    }
+  return 0;
+}
+
+
+/* Remove leading blanks from a string, padding at end.  The src and dest
+   should not overlap.  */
+
+void
+adjustl (char *dest, GFC_INTEGER_4 len, const char *src)
+{
+  int i;
+
+  i = 0;
+  while (i<len && src[i] == ' ')
+    i++;
+
+  if (i < len)
+    memcpy (dest, &src[i], len - i);
+  if (i > 0)
+    memset (&dest[len - i], ' ', i);
+}
+
+
+/* Remove trailing blanks from a string.  */
+
+void
+adjustr (char *dest, GFC_INTEGER_4 len, const char *src)
+{
+  int i;
+
+  i = len;
+  while (i > 0 && src[i - 1] == ' ')
+    i++;
+
+  if (i < len)
+    memcpy (&dest[len - i], &src, i);
+  if (i < len)
+    memset (dest, ' ', len - i);
+}
+
+
+/* Scan a string for any one of the characters in a set of characters.  */
+
+GFC_INTEGER_4
+string_scan (GFC_INTEGER_4 slen, const char * str, GFC_INTEGER_4 setlen,
+             const char * set, GFC_LOGICAL_4 back)
+{
+  int start;
+  int last;
+  int i;
+  int delta;
+
+  if (slen == 0 || setlen == 0)
+    return 0;
+
+  if (back)
+    {
+      last =  0;
+      start = slen - 1;
+      delta = -1;
+    }
+  else
+    {
+      last = slen - 1;
+      start = 0;
+      delta = 1;
+    }
+
+  i = 0;
+  for (; start != last; start += delta)
+    {
+      for (i = 0; i < setlen; i++)
+        {
+          if (str[start] == set[i])
+            return (start + 1);
+        }
+    }
+
+  return 0;
+}
+
+
+/* Verify that a set of characters contains all the characters in a
+   string by indentifying the position of the first character in a
+   characters that dose not appear in a given set of characters.  */
+
+GFC_INTEGER_4
+string_verify (GFC_INTEGER_4 slen, const char * str, GFC_INTEGER_4 setlen,
+               const char * set, GFC_LOGICAL_4 back)
+{
+  int start;
+  int last;
+  int i;
+  int delta;
+
+  if (slen == 0)
+    return 0;
+
+  if (back)
+    {
+      last =  0;
+      start = slen - 1;
+      delta = -1;
+    }
+  else
+    {
+      last = slen - 1;
+      start = 0;
+      delta = 1;
+    }
+  i = 0;
+  for (; start != last; start += delta)
+    {
+      for (i = 0; i < setlen; i++)
+        {
+          if (str[start] == set[i])
+            break;
+        }
+      if (i == setlen)
+        return (start + 1);
+    }
+
+  return 0;
+}
+
+
+/* Concatenate several copies of a string.  */
+
+void
+string_repeat (char * dest, GFC_INTEGER_4 slen, 
+               const char * src, GFC_INTEGER_4 ncopies)
+{
+  int i;
+
+  /* See if ncopies is valid.  */
+  if (ncopies < 0)
+    {
+      /* The error is already reported.  */
+      runtime_error ("Augument NCOPIES is negative.");
+    }
+
+  /* Copy characters.  */
+  for (i = 0; i < ncopies; i++) 
+    {
+      memmove (dest + (i * slen), src, slen);
+    }
+}
+
diff --git a/libgfortran/intrinsics/transpose_generic.c b/libgfortran/intrinsics/transpose_generic.c
new file mode 100644
index 00000000000..d72ae5a4b81
--- /dev/null
+++ b/libgfortran/intrinsics/transpose_generic.c
@@ -0,0 +1,74 @@
+/* Implementation of the TRANSPOSE intrinsic
+   Copyright 2003 Free Software Foundation, Inc.
+   Contributed by Tobias Schlüter
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfor is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Ligbfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "config.h"
+#include <stdlib.h>
+#include <string.h>
+#include <assert.h>
+#include "libgfortran.h"
+
+void
+__transpose (gfc_array_char * ret, gfc_array_char * source)
+{
+  /* r.* indicates the return array.  */
+  index_type rxstride, rystride;
+  char *rptr;
+  /* s.* indicates the source array.  */
+  index_type sxstride, systride;
+  const char *sptr;
+
+  index_type xcount, ycount;
+  index_type x, y;
+  index_type size;
+
+  assert (GFC_DESCRIPTOR_RANK (source) == 2
+          && GFC_DESCRIPTOR_RANK (ret) == 2);
+
+  size = GFC_DESCRIPTOR_SIZE (source);
+  sxstride = source->dim[0].stride * size;
+  if (sxstride == 0)
+    sxstride = size;
+  systride = source->dim[1].stride * size;
+  xcount = source->dim[0].ubound + 1 - source->dim[0].lbound;
+  ycount = source->dim[1].ubound + 1 - source->dim[1].lbound;
+
+  rxstride = ret->dim[0].stride * size;
+  if (rxstride == 0)
+    rxstride = size;
+  rystride = ret->dim[1].stride * size;
+
+  rptr = ret->data;
+  sptr = source->data;
+
+  for (y = 0; y < ycount; y++)
+    {
+      for (x = 0; x < xcount; x++)
+        {
+          memcpy (rptr, sptr, size);
+
+          sptr += sxstride;
+          rptr += rystride;
+        }
+      sptr += systride - (sxstride * xcount);
+      rptr += rxstride - (rystride * xcount);
+    }
+}
+
diff --git a/libgfortran/intrinsics/unpack_generic.c b/libgfortran/intrinsics/unpack_generic.c
new file mode 100644
index 00000000000..3d02a3e6060
--- /dev/null
+++ b/libgfortran/intrinsics/unpack_generic.c
@@ -0,0 +1,154 @@
+/* Generic implementation of the RESHAPE intrinsic
+   Copyright 2002 Free Software Foundation, Inc.
+   Contributed by Paul Brook <paul@nowt.org>
+
+This file is part of the GNU Fortran 95 runtime library (libgfor).
+
+Libgfor is free software; you can redistribute it and/or
+modify it under the terms of the GNU Lesser General Public
+License as published by the Free Software Foundation; either
+version 2.1 of the License, or (at your option) any later version.
+
+Ligbfor is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public
+License along with libgfor; see the file COPYING.LIB.  If not,
+write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA.  */
+
+#include "config.h"
+#include <stdlib.h>
+#include <assert.h>
+#include <string.h>
+#include "libgfortran.h"
+
+void
+__unpack1 (const gfc_array_char * ret, const gfc_array_char * vector,
+    const gfc_array_l4 * mask, const gfc_array_char * field)
+{
+  /* r.* indicates the return array.  */
+  index_type rstride[GFC_MAX_DIMENSIONS];
+  index_type rstride0;
+  char *rptr;
+  /* v.* indicates the vector array.  */
+  index_type vstride0;
+  char *vptr;
+  /* f.* indicates the field array.  */
+  index_type fstride[GFC_MAX_DIMENSIONS];
+  index_type fstride0;
+  const char *fptr;
+  /* m.* indicates the mask array.  */
+  index_type mstride[GFC_MAX_DIMENSIONS];
+  index_type mstride0;
+  const GFC_LOGICAL_4 *mptr;
+
+  index_type count[GFC_MAX_DIMENSIONS];
+  index_type extent[GFC_MAX_DIMENSIONS];
+  index_type n;
+  index_type dim;
+  index_type size;
+  index_type fsize;
+
+  size = GFC_DESCRIPTOR_SIZE (ret);
+  /* A field element size of 0 actually means this is a scalar.  */
+  fsize = GFC_DESCRIPTOR_SIZE (field);
+  dim = GFC_DESCRIPTOR_RANK (ret);
+  for (n = 0; n < dim; n++)
+    {
+      count[n] = 0;
+      extent[n] = ret->dim[n].ubound + 1 - ret->dim[n].lbound;
+      rstride[n] = ret->dim[n].stride * size;
+      fstride[n] = field->dim[n].stride * fsize;
+      mstride[n] = mask->dim[n].stride;
+    }
+  if (rstride[0] == 0)
+    rstride[0] = size;
+  if (fstride[0] == 0)
+    fstride[0] = fsize;
+  if (mstride[0] == 0)
+    mstride[0] = 1;
+
+  vstride0 = vector->dim[0].stride * size;
+  if (vstride0 == 0)
+    vstride0 = size;
+  rstride0 = rstride[0];
+  fstride0 = fstride[0];
+  mstride0 = mstride[0];
+  rptr = ret->data;
+  fptr = field->data;
+  mptr = mask->data;
+  vptr = vector->data;
+
+
+  /* Use the same loop for both logical types. */
+  if (GFC_DESCRIPTOR_SIZE (mask) != 4)
+    {
+      if (GFC_DESCRIPTOR_SIZE (mask) != 8)
+        runtime_error ("Funny sized logical array");
+      for (n = 0; n < dim; n++)
+        mstride[n] <<= 1;
+      mstride0 <<= 1;
+      mptr = GFOR_POINTER_L8_TO_L4 (mptr);
+    }
+
+  while (rptr)
+    {
+      if (*mptr)
+        {
+          /* From vector.  */
+          memcpy (rptr, vptr, size);
+          vptr += vstride0;
+        }
+      else
+        {
+          /* From field.  */
+          memcpy (rptr, fptr, size);
+        }
+      /* Advance to the next element.  */
+      rptr += rstride0;
+      fptr += fstride0;
+      mptr += mstride0;
+      count[0]++;
+      n = 0;
+      while (count[n] == extent[n])
+        {
+          /* When we get to the end of a dimension, reset it and increment
+             the next dimension.  */
+          count[n] = 0;
+          /* We could precalculate these products, but this is a less
+             frequently used path so proabably not worth it.  */
+          rptr -= rstride[n] * extent[n];
+          fptr -= fstride[n] * extent[n];
+          mptr -= mstride[n] * extent[n];
+          n++;
+          if (n >= dim)
+            {
+              /* Break out of the loop.  */
+              rptr = NULL;
+              break;
+            }
+          else
+            {
+              count[n]++;
+              rptr += rstride[n];
+              fptr += fstride[n];
+              mptr += mstride[n];
+            }
+        }
+    }
+}
+
+void
+__unpack0 (const gfc_array_char * ret, const gfc_array_char * vector,
+    const gfc_array_l4 * mask, char * field)
+{
+  gfc_array_char tmp;
+
+  tmp.dtype = 0;
+  tmp.data = field;
+  __unpack1 (ret, vector, mask, &tmp);
+}
+