[Bug fortran/36158] New: Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[Bug fortran/36158] New: Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

Fortran 2008 (PR 33197) supports the BESSEL_* functions. However, for BESSEL_*N
there are two versions:

A) An elemental function BESSEL_YN(N, X) and BESSEL_JN(N,X)
B) An transformational function BESSEL_YN(N1,N2,X) and BESSEL_JN(N1,N2,Y)

(B) Is currently not supported:
The result of BESSEL_JN (N1, N2, X) is a rank-one array with extent
MAX(N2−N1+1,0). Element i of the result value of BESSEL_JN (N1, N2, X) is
a processor-dependent approximation to the Bessel function of the first kind
and order N1+i − 1 of X.

Note that one needs to simplify also the transformational version as one can
also use it in initialization expressions.

See also PR 36096 and 36117 and
http://groups.google.com/group/comp.lang.fortran/browse_thread/thread/e9add97708681397


-- 
           Summary: Transformational function BESSEL_YN(n1,n2,x) and
                    BESSEL_JN missing
           Product: gcc
           Version: 4.4.0
            Status: UNCONFIRMED
          Severity: normal
          Priority: P3
         Component: fortran
        AssignedTo: unassigned at gcc dot gnu dot org
        ReportedBy: burnus at gcc dot gnu dot org


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[fxcoudert at gcc dot gnu dot org]

-- 

fxcoudert at gcc dot gnu dot org changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
             Status|UNCONFIRMED                 |NEW
     Ever Confirmed|0                           |1
   Last reconfirmed|0000-00-00 00:00:00         |2009-03-28 15:46:40
               date|                            |


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[fxcoudert at gcc dot gnu dot org]

-- 

fxcoudert at gcc dot gnu dot org changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
           Severity|normal                      |enhancement


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

------- Comment #1 from burnus at gcc dot gnu dot org  2010-08-16 19:33 -------
Created an attachment (id=21496)
 --> (http://gcc.gnu.org/bugzilla/attachment.cgi?id=21496&action=view)
First patch

First implementation with compile-time resolution

Todo:
- Handle run-time solution - I do not yet know how to handle it best.
  One solution is to use effectively  [ (bessel_jn(i,x), i=n1, n2) ]
  but the question is where; maybe in simplify, but then one needs to
  generate a temporary symbol which does not clash with others.
  Cf. array.c's match_array_list for an ac-implied-do
- Test cases
- Documentation (cf. also PR 36096)


-- 


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

------- Comment #2 from burnus at gcc dot gnu dot org  2010-08-19 07:28 -------
Subject: Bug 36158

Author: burnus
Date: Thu Aug 19 07:28:17 2010
New Revision: 163364

URL: http://gcc.gnu.org/viewcvs?root=gcc&view=rev&rev=163364
Log:
2010-08-19  Tobias Burnus  <burnus@net-b.de>

        PR fortran/36158
        PR fortran/33197
        * check.c (gfc_check_bessel_n2): New function.
        * gfortran.h (gfc_isym_id): Add GFC_ISYM_JN2 and GFC_ISYM_YN2.
        * intrinsic.c (add_functions): Add transformational version
        of the Bessel_jn/yn intrinsics.
        * intrinsic.h (gfc_check_bessel_n2,gfc_simplify_bessel_jn2,
        gfc_simplify_bessel_yn2): New prototypes.
        * intrinsic.texi (Bessel_jn, Bessel_yn): Document
        transformational variant.
        * simplify.c (gfc_simplify_bessel_jn, gfc_simplify_bessel_yn):
        Check for negative order.
        (gfc_simplify_bessel_n2,gfc_simplify_bessel_jn2,
        gfc_simplify_bessel_yn2): New functions.

2010-08-19  Tobias Burnus  <burnus@net-b.de>

        PR fortran/36158
        PR fortran/33197
        * gfortran.dg/bessel_3.f90: New.
        * gfortran.dg/bessel_4.f90: New.
        * gfortran.dg/bessel_5.f90: New.


Added:
    trunk/gcc/testsuite/gfortran.dg/bessel_3.f90
    trunk/gcc/testsuite/gfortran.dg/bessel_4.f90
    trunk/gcc/testsuite/gfortran.dg/bessel_5.f90
Modified:
    trunk/gcc/fortran/ChangeLog
    trunk/gcc/fortran/check.c
    trunk/gcc/fortran/gfortran.h
    trunk/gcc/fortran/intrinsic.c
    trunk/gcc/fortran/intrinsic.h
    trunk/gcc/fortran/intrinsic.texi
    trunk/gcc/fortran/simplify.c
    trunk/gcc/testsuite/ChangeLog


-- 


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

------- Comment #3 from burnus at gcc dot gnu dot org  2010-08-19 07:37 -------
The patch implemented the compile-time version using a recurrence algorithm.

TODO: Add a run-time implementation (using a recurrence algorithm).


-- 


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

------- Comment #4 from burnus at gcc dot gnu dot org  2010-08-19 08:21 -------
(In reply to comment #3)
> TODO: Add a run-time implementation (using a recurrence algorithm).

Remember that all kinds of integers are allowed for N1 and N2, but that POSIX's
jn/yn only allows "int" - thus, on can add a cast to c_int in iresolve.c.


-- 


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

------- Comment #5 from burnus at gcc dot gnu dot org  2010-08-19 23:08 -------
Created an attachment (id=21526)
 --> (http://gcc.gnu.org/bugzilla/attachment.cgi?id=21526&action=view)
Run-time implementation

Implementation works, but I would like to pass the arguments by VALUE not by
reference. Additionally, a test case is missing.


-- 

burnus at gcc dot gnu dot org changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
  Attachment #21496|0                           |1
        is obsolete|                            |


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[kargl at gcc dot gnu dot org]

------- Comment #6 from kargl at gcc dot gnu dot org  2010-08-19 23:28 -------
(In reply to comment #5)
> Created an attachment (id=21526)
 --> (http://gcc.gnu.org/bugzilla/attachment.cgi?id=21526&action=view) [edit]
> Run-time implementation
> 
> Implementation works, but I would like to pass the arguments by VALUE not by
> reference. Additionally, a test case is missing.
> 

I read the patch.  It looks correct to me.  This
is OK to commit ever with the by reference semantics.


-- 


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

------- Comment #7 from burnus at gcc dot gnu dot org  2010-08-20 08:43 -------
Comparison of the recurrence version with calling the library every time
(for the run-time library, I expect a similar result for MPFR/compile time
version):
http://gcc.gnu.org/ml/fortran/2010-08/msg00252.html


-- 


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

------- Comment #8 from burnus at gcc dot gnu dot org  2010-08-20 17:23 -------
Created an attachment (id=21532)
 --> (http://gcc.gnu.org/bugzilla/attachment.cgi?id=21532&action=view)
Run-time implementation

Possibly final patch, though not yet regtested.


-- 

burnus at gcc dot gnu dot org changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
  Attachment #21526|0                           |1
        is obsolete|                            |


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

------- Comment #9 from burnus at gcc dot gnu dot org  2010-08-20 17:50 -------
(In reply to comment #8)
> Created an attachment (id=21532)
 --> (http://gcc.gnu.org/bugzilla/attachment.cgi?id=21532&action=view) [edit]
> Possibly final patch, though not yet regtested.

Fails for:
  integer :: x(2)
  x = [1,2]
  print *, RESHAPE([x], [1,2])
  end

Probably due to trans-intrinsics.c change of copying the formal arg list - even
though it looks OK in the debugger. I could not yet trace down where exactly it
fails. GDB/valgrind show the line:

 Invalid read of size 4
    at 0x57A7F9: gfc_conv_procedure_call (trans-expr.c:2755)
    by 0x58AF4E: conv_generic_with_optional_char_arg.constprop.32
(trans-intrinsic.c:3326)
    by 0x58BBDA: gfc_conv_intrinsic_function (trans-intrinsic.c:4917)

but trans-expr.c:2755 is the line where gfc_conv_procedure_call is defined,
which does not help me with debugging :-(


-- 


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

------- Comment #10 from burnus at gcc dot gnu dot org  2010-08-21 10:13 -------
Subject: Bug 36158

Author: burnus
Date: Sat Aug 21 10:12:53 2010
New Revision: 163440

URL: http://gcc.gnu.org/viewcvs?root=gcc&view=rev&rev=163440
Log:
2010-08-21  Tobias Burnus  <burnus@net-b.de>

        PR fortran/36158
        PR fortran/33197
        * intrinsic.c (add_sym): Init value attribute.
        (set_attr_value): New function.
        (add_functions) Use it and add JN/YN resolvers.
        * symbol.c (gfc_copy_formal_args_intr): Copy value attr.
        * intrinsic.h (gfc_resolve_bessel_n2): New prototype.
        * gfortran.h (gfc_intrinsic_arg): Add value attribute.
        * iresolve.c (gfc_resolve_bessel_n2): New function.
        * trans-intrinsic.c (gfc_get_symbol_for_expr): Create
        formal arg list.
        (gfc_conv_intrinsic_function,gfc_is_intrinsic_libcall):
        Add GFC_ISYM_JN2/GFC_ISYM_YN2 as case value.
        * simplify.c (): For YN set to -INF if previous values
        was -INF.
        * trans-expr.c (gfc_conv_procedure_call): Don't crash
        if sym->as is NULL.
        * iresolve.c (gfc_resolve_extends_type_of): Set the
        type of the dummy argument to the one of the actual.

2010-08-21  Tobias Burnus  <burnus@net-b.de>

        PR fortran/36158
        PR fortran/33197
        * m4/bessel.m4: Implement bessel_jn and bessel_yn.
        * gfortran.map: Add the generated bessel_jn_r{4,8,10,16}
        and bessel_yn_r{4,8,10,16}.
        * Makefile.am: Add bessel.m4.
        * Makefile.in: Regenerated.
        * generated/bessel_r4.c: Generated.
        * generated/bessel_r16.c: Generated.
        * generated/bessel_r8.c: Generated.
        * generated/bessel_r10.c: Generated.

2010-08-21  Tobias Burnus  <burnus@net-b.de>

        PR fortran/36158
        PR fortran/33197
        * gfortran.dg/bessel_6.f90: New.
        * gfortran.dg/bessel_7.f90: New.


Added:
    trunk/gcc/testsuite/gfortran.dg/bessel_6.f90
    trunk/gcc/testsuite/gfortran.dg/bessel_7.f90
    trunk/libgfortran/generated/bessel_r10.c
    trunk/libgfortran/generated/bessel_r16.c
    trunk/libgfortran/generated/bessel_r4.c
    trunk/libgfortran/generated/bessel_r8.c
    trunk/libgfortran/m4/bessel.m4
Modified:
    trunk/gcc/fortran/ChangeLog
    trunk/gcc/fortran/gfortran.h
    trunk/gcc/fortran/intrinsic.c
    trunk/gcc/fortran/intrinsic.h
    trunk/gcc/fortran/iresolve.c
    trunk/gcc/fortran/simplify.c
    trunk/gcc/fortran/symbol.c
    trunk/gcc/fortran/trans-expr.c
    trunk/gcc/fortran/trans-intrinsic.c
    trunk/gcc/testsuite/ChangeLog
    trunk/libgfortran/ChangeLog
    trunk/libgfortran/Makefile.am
    trunk/libgfortran/Makefile.in
    trunk/libgfortran/gfortran.map


-- 


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

------- Comment #11 from burnus at gcc dot gnu dot org  2010-08-21 10:28 -------
FIXED on the 4.6 trunk.


-- 

burnus at gcc dot gnu dot org changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
             Status|NEW                         |RESOLVED
         Resolution|                            |FIXED


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[dominiq at lps dot ens dot fr]

------- Comment #12 from dominiq at lps dot ens dot fr  2010-08-21 13:11 -------
On x86_64-apple-darwin10.4 I need the following adjustments:

--- /opt/gcc/work/gcc/testsuite/gfortran.dg/bessel_6.f90        2010-08-21
12:30:29.000000000 +0200
+++ bessel_6_db.f90     2010-08-21 15:05:20.000000000 +0200
@@ -8,7 +8,7 @@
 implicit none
 real,parameter :: values(*) = [0.0, 0.5, 1.0, 0.9,
1.8,2.0,3.0,4.0,4.25,8.0,34.53, 475.78] 
 real,parameter :: myeps(size(values)) = epsilon(0.0) &
-                  * [2, 7, 5, 6, 9, 12, 12, 7, 7, 8, 75, 6 ]
+                  * [2, 7, 5, 6, 9, 12, 12, 7, 7, 8, 120, 15 ]
 ! The following is sufficient for me - the values above are a bit
 ! more tolerant
 !                  * [0, 5, 3, 4, 6, 7, 7, 5, 5, 6, 66, 4 ]

--- /opt/gcc/work/gcc/testsuite/gfortran.dg/bessel_7.f90        2010-08-21
12:30:30.000000000 +0200
+++ bessel_7_db.f90     2010-08-21 15:02:49.000000000 +0200
@@ -8,7 +8,7 @@
 implicit none
 real,parameter :: values(*) = [0.0, 0.5, 1.0, 0.9,
1.8,2.0,3.0,4.0,4.25,8.0,34.53, 475.78] 
 real,parameter :: myeps(size(values)) = epsilon(0.0) &
-                  * [2, 2, 2, 5, 5, 2, 12, 2, 4, 3, 30, 130 ]
+                  * [2, 3, 2, 5, 7, 2, 12, 3, 7, 3, 30, 144 ]
 ! The following is sufficient for me - the values above are a bit
 ! more tolerant
 !                  * [0, 0, 0, 3, 3, 0, 9, 0, 2, 1, 22, 130 ]
@@ -19,7 +19,7 @@ real :: rec(0:Nmax), lib(0:Nmax)
 integer :: i

 do i = 1, ubound(values,dim=1)
-  call compare(values(i), myeps(i), nit(i), 3*epsilon(0.0))
+  call compare(values(i), myeps(i), nit(i), myeps(i))
 end do

 contains


-- 


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

------- Comment #13 from burnus at gcc dot gnu dot org  2010-08-22 08:34 -------
Subject: Bug 36158

Author: burnus
Date: Sun Aug 22 08:34:14 2010
New Revision: 163454

URL: http://gcc.gnu.org/viewcvs?root=gcc&view=rev&rev=163454
Log:
2010-08-22  Tobias Burnus  <burnus@net-b.de>
            Dominique d'Humieres <dominiq@lps.ens.fr>

        PR fortran/45367
        PR fortran/36158
        * gfortran.dg/bessel_6.f90: Increase numeric tolerence.
        * gfortran.dg/bessel_7.f90: Increase numeric tolerence.


Modified:
    trunk/gcc/testsuite/ChangeLog
    trunk/gcc/testsuite/gfortran.dg/bessel_7.f90


-- 


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[burnus at gcc dot gnu dot org]

------- Comment #14 from burnus at gcc dot gnu dot org  2010-08-22 08:40 -------
Subject: Bug 36158

Author: burnus
Date: Sun Aug 22 08:39:48 2010
New Revision: 163455

URL: http://gcc.gnu.org/viewcvs?root=gcc&view=rev&rev=163455
Log:
2010-08-22  Tobias Burnus  <burnus@net-b.de>

        PR fortran/36158
        * gfortran.dg/bessel_7.f90: Disable accidently enabled debug
        output.


Modified:
    trunk/gcc/testsuite/ChangeLog
    trunk/gcc/testsuite/gfortran.dg/bessel_7.f90


-- 


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[ibolton at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

--- Comment #17 from Ian Bolton <ibolton at gcc dot gnu.org> 2011-02-07 15:41:14 UTC ---
(In reply to comment #16)
> (In reply to comment #15)
> 
> > 
> > Is this going to be backported to 4.5?  It is required to make CSHIFT function
> > correctly.
> 
> It may be prudent to actually describe the problem
> with cshift and give a short example.  The above
> is of little value.

Sorry about that.  You are right - I should have given more information.  Here
goes:

Currently, the code for the gfc_get_symbol_for_expr() function on the 4.5
branch has a comment that says:

"/* TODO: proper argument lists for external intrinsics.  */"

CSHIFT is one such fortran intrinsic, which will have an incorrect formal
argument list due to this TODO not being done.  Specifically, we can end up
with a _gfortran_cshift0_4 call with no formal arguments.

Here is a reduced testcase (adapted from an existing testcase) that shows such
an issue with _gfortran_cshift0_4:

 ! Program to test the cshift intrinsic
 program intrinsic_cshift
    integer, dimension(3, 2) :: a

    ! Scalar shift
    a = reshape ((/1, 2, 3, 4, 5, 6/), (/3, 2/))
    a = cshift (a, 1, 1)
    if (any (a .ne. reshape ((/2, 3, 1, 5, 6, 4/), (/3, 2/)))) &
       call abort
 end program

There is no information within the existing tree dumps to show things have gone
wrong, but debugging with gdb and inserting a breakpoint within the
gfc_get_symbol_for_expr function will show that
expr->value.function.isym.formal has three items in it, whereas the new sym
that we create is empty.  The patch that fixed PR36158 correctly copies the
formal arguments from isym to sym and we then get the correct behaviour.

I hope that information is sufficient.

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[kargl at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

kargl at gcc dot gnu.org changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
                 CC|                            |kargl at gcc dot gnu.org

--- Comment #16 from kargl at gcc dot gnu.org 2011-02-07 15:21:11 UTC ---
(In reply to comment #15)

> 
> Is this going to be backported to 4.5?  It is required to make CSHIFT function
> correctly.

It may be prudent to actually describe the problem
with cshift and give a short example.  The above
is of little value.

[Bug fortran/36158] Transformational function BESSEL_YN(n1,n2,x) and BESSEL_JN missing

[ibolton at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36158

Ian Bolton <ibolton at gcc dot gnu.org> changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
                 CC|                            |ibolton at gcc dot gnu.org

--- Comment #15 from Ian Bolton <ibolton at gcc dot gnu.org> 2011-02-07 13:47:35 UTC ---
(In reply to comment #10)
> Subject: Bug 36158
> 
> Author: burnus
> Date: Sat Aug 21 10:12:53 2010
> New Revision: 163440
> 
> URL: http://gcc.gnu.org/viewcvs?root=gcc&view=rev&rev=163440
> Log:
> 2010-08-21  Tobias Burnus  <burnus@net-b.de>
> 
>         PR fortran/36158
>         PR fortran/33197
>         * intrinsic.c (add_sym): Init value attribute.
>         (set_attr_value): New function.
>         (add_functions) Use it and add JN/YN resolvers.
>         * symbol.c (gfc_copy_formal_args_intr): Copy value attr.
>         * intrinsic.h (gfc_resolve_bessel_n2): New prototype.
>         * gfortran.h (gfc_intrinsic_arg): Add value attribute.
>         * iresolve.c (gfc_resolve_bessel_n2): New function.
>         * trans-intrinsic.c (gfc_get_symbol_for_expr): Create
>         formal arg list.
>         (gfc_conv_intrinsic_function,gfc_is_intrinsic_libcall):
>         Add GFC_ISYM_JN2/GFC_ISYM_YN2 as case value.
>         * simplify.c (): For YN set to -INF if previous values
>         was -INF.
>         * trans-expr.c (gfc_conv_procedure_call): Don't crash
>         if sym->as is NULL.
>         * iresolve.c (gfc_resolve_extends_type_of): Set the
>         type of the dummy argument to the one of the actual.
> 
> 2010-08-21  Tobias Burnus  <burnus@net-b.de>
> 
>         PR fortran/36158
>         PR fortran/33197
>         * m4/bessel.m4: Implement bessel_jn and bessel_yn.
>         * gfortran.map: Add the generated bessel_jn_r{4,8,10,16}
>         and bessel_yn_r{4,8,10,16}.
>         * Makefile.am: Add bessel.m4.
>         * Makefile.in: Regenerated.
>         * generated/bessel_r4.c: Generated.
>         * generated/bessel_r16.c: Generated.
>         * generated/bessel_r8.c: Generated.
>         * generated/bessel_r10.c: Generated.
> 
> 2010-08-21  Tobias Burnus  <burnus@net-b.de>
> 
>         PR fortran/36158
>         PR fortran/33197
>         * gfortran.dg/bessel_6.f90: New.
>         * gfortran.dg/bessel_7.f90: New.
> 


Is this going to be backported to 4.5?  It is required to make CSHIFT function
correctly.