In the following algorithm I set up a do loop in such a way to perform some operations for each value of the variable domega. The problem of my code is that it misses the last value of domega, that should be equal to tau2.

Can you tell me how to solve this issue? In matlab I would write an array of domega like: domega = first_value:step:last_value and then I would perform a for loop on the values of domega array.

Can you show how to modify my code to get a do loop on domega, by starting from first value up to last value with the indicated step size?

  domega  = tau1                       ! first value 
   
  if(tau2 < 0.d0) tau2 = tau2   twopi  ! tau2 is the last value
      
  delta    =  4.617207817725326D-003  ! Interval step size
  n_domega = 1   ((tau2 - tau1)/delta)     ! numbers of values within interval


  do i_om = 1,n_domega    ! Loop 
      
     
 !  ( COMPUTATIONS...)
   
   
      
      domega = tau1   i_om*delta
  enddo  ! end loop on Delta_omega 

PROBLEM SOLVED! After reading your comments and some web searches, I managed to solve the problem by using the following linspace subroutine:

!***********************************************************************
    subroutine linspace_real(from, to, array)
    !***********************************************************************
    ! Generates evenly spaced numbers from `from` to `to` (inclusive).
    !
    ! Inputs:
    ! -------
    !
    ! from, to : the lower and upper boundaries of the numbers to generate
    !
    ! Outputs:
    ! -------
    !
    ! array : Array of evenly spaced numbers
    !***********************************************************************
        implicit none
    
        !Arguments
        real(pr), intent(in)  :: from, to
        real(pr), intent(out) :: array(:)
        
        real(pr) :: range
        integer :: n, i
    
        n = size(array)
        range = to - from

        if (n == 0) return

        if (n == 1) then
            array(1) = from
            return
        end if


        do i=1, n
            array(i) = from   range * (i - 1) / (n - 1)
        end do
        
        return
        
    end subroutine linspace_real
    !*********************************************************************** 

In the main program I wrote:

!Declarations
real(pr), dimension(:), allocatable  :: d_omega_vec 

 
! Scan interval of Delta_omega (we know the interval step, so we have to compute the number of values within the interval)
domega  = tau1
                     
if(tau2 < 0.d0) tau2 = tau2   twopi
                      
delta    =  4.617207817725326D-003  ! Interval step size
n_domega = 1   ((tau2 - tau1)/delta)     ! numbers of values within interval
              
n_domega_vec = n_domega   1 
                    
allocate(d_omega_vec(n_domega_vec))
call linspace_real(tau1, tau2, d_omega_vec)

So, instead to increment domega during the loop, I created an array of domega from tau1 to tau2 with number of elements properly computed by exploiting my delta.

CodePudding user response：

I'm not sure if delta = 0.0026 (as in the code) or delta= 4.617207817725326D-003 (as in the comment), but here is a way to give an initial value and loop adding a step each time:

program test
  double precision, parameter :: tau1 = 0.392464539129341, tau2 = 5.89072076805024, delta = 0.0026 
  double precision :: domega

  domega  = tau1                       ! first value 
  do
    if(domega > tau2) exit
    write(*,"('domega = ',es16.9)") domega
    domega = domega   delta
  enddo  ! end loop on Delta_omega 

  write(*,"('next domega would be = ',es16.9)") domega   delta
end program test

End of output:

(...)
domega =  5.883664563E 00
domega =  5.886264563E 00
domega =  5.888864563E 00
next domega would be =  5.894064563E 00

The last line is printed outside of the loop, and is already above tau2 = 5.89072076805024.

CodePudding user response：

From what I can tell what you are looking for is basically Python's numpy.linspace. You can write effectively the same in Fortran with an implied do-loop:

array = start_val   val_spacing * [(i, i=0,N-1)]

Sample

program main
   use :: iso_fortran_env, only: dp=>real64
   implicit none

   real(dp), allocatable :: arr(:)

   arr =  linspace(1.0_dp, 2.0_dp, 11)
   print*, arr

   contains

   function linspace(a,b,N) result(lin_arr)
        !! Equivalent to Numpy's numpy.linspace

        real(dp), intent(in) :: a,b
        integer, intent(in) :: N
        real(dp), allocatable :: lin_arr(:)

        ! Local
        integer :: i
        real(dp) :: dx

        dx = (b-a) / (N-1)
        lin_arr = a   dx * [(i, i=0, N-1)]
   end function

end program main