1

I am having trouble in implementing an approach to call Newton's method in a Fortran program. So I want to use Newton's method to solve an equation following the link

However, my program is slightly different with the example above. In my case, the equation requires some additional information which are produced during runtime.

subroutine solve(f, fp, x0, x, iters, debug)

which means the f is calculated not only based on x, but also a few other variables (but x is the unknown).

I have a solution, which only works for a simple case: I used a module to include Newton's solver. I also defined a derived data type to hold all the arguments inside the module. It works good now.

My question is: I need to call the Newton's method many times, and each time the arguments are different. How should I design the structure of the modules? Or should I use a different solution?

I provided a simple example below:

module solver
  type argu
    integer :: m
  end type argu
  type(argu):: aArgu_test  !should I put here?
  contains
    subroutine solve(f, fp, x0, x, iters, debug)
       ...
       !m is used inside here
    end subroutine solve
    subroutine set_parameter(m_in)
       aArgu%m = m_in
    end subroutine set_parameter()

end module solver

And the calling module is:

!only one set of argument, but used many times
module A
   use module solver

   do i = 1, 4, 1
     set_parameter(i) 
     !call Newtow method
     ...
   enddo
end module A

!can I use an array for argu type if possible?
module B
   use module solver
   type(argu), dimension(:), allocable :: aArgu ! or should I put here or inside solver module?

end module B

My understanding is that if I put the argu object inside the solver module, then all solver calling will use the same arguments (I can still save all of them inside module A using the above method). In that case, I have to update the arguments during each for loop?

Because the program runs using MPI/OpenMP, I want to make sure there is no overwritten among threads. Thank you.

Vladimir F Героям слава
  • 57,977
  • 4
  • 76
  • 119
Chang
  • 396
  • 4
  • 17
  • When you say "call the Newton's method many times", do you mean: find the target `x` for a given set of parameters, and then find a new target `x` for a different set of parameters? – Ross Sep 06 '18 at 04:27
  • Potential duplicates https://stackoverflow.com/questions/24127313/passing-external-function-of-multiple-variables-as-a-function-of-one-variable-in?noredirect=1&lq=1 https://stackoverflow.com/questions/26297170/fortran-minimization-of-a-function-with-additional-arguments https://stackoverflow.com/questions/38324934/function-with-more-arguments-and-integration?noredirect=1&lq=1 https://stackoverflow.com/questions/37714406/implementing-anonymous-functions-in-fortran?noredirect=1&lq=1 – Vladimir F Героям слава Sep 06 '18 at 06:21
  • Please keep your questions to a single question only. But it seems you actually answered your qestion 1 yourself... But we need more details for your question 2. Please show us how does the code actually work, simply put - show more code. – Vladimir F Героям слава Sep 06 '18 at 06:24
  • So not closing yet, because I assume you are actually asking for Q2 and Q1 is solved. But we need more info. – Vladimir F Героям слава Sep 06 '18 at 06:26
  • Well, I want to use the Newton's method in my program, the first step/question is a test with additional arguments. My final program requires more than one sets of arguments. So @Ross is right. I will find x for multiple sets of arguments. – Chang Sep 06 '18 at 16:53
  • Have you looked at the answers linked by @VladimirF? It looks to me like they answer your question. – Ross Sep 06 '18 at 16:57
  • @VladimirF provided some useful information, which are very close to my question 1. I have edited the part for question 2. – Chang Sep 06 '18 at 17:31

2 Answers2

3

There is a common pattern in modern Fortran for the problem you are facing (partial function application). Unlike other languages, Fortran doesn't have function closures, so making a lexical scope for a function is a little "convoluted" and kind of limited.

You should really consider revisiting all the links @VladmirF shared on the comment, most of them apply straightforwardly to your case. I will give you an example of a solution.

This is a solution without using a wrapper type. I will use a feature included in Fortran 2008 standard: passing an internal procedure as an argument. It is compatible with the latest gfortran, Intel and many others. If you can't access a compiler with this feature or if you prefer a solution with a derived type, you can refer to this answer.

module without_custom_type
  use, intrinsic :: iso_fortran_env, only: r8 => real64
  use :: solver

contains
  subroutine solve_quad(a, b, c, x0, x, iters, debug)
    integer, intent(in) :: a, b, c
    real(r8), intent(in) :: x0
    real(r8), intent(out) :: x
    integer, intent(out) :: iters
    logical, intent(in) :: debug

    call solve(f, fp, x0, x, iters, debug)

  contains
    real(r8) function f(x)
      real(r8),intent(in) :: x
      f = a * x * x + b * x + c
    end

    real(r8) function fp(x)
      real(r8),intent(in) :: x
      fp = 2 * a * x + b
    end
  end
end

The rationale of this code is: as f and fp lay inside of the solve_quad procedure, they have access to the arguments a, b and c by host association, without touching those function's signatures. The resulting effect is like changing the arity of the function.

Testing it with gfortran 8.0 and the solver implementation from the link you shared, I got this:

program test
  use, intrinsic :: iso_fortran_env, only: r8 => real64
  use :: without_custom_type
  implicit none

  real(r8) :: x, x0
  integer :: iters
  integer :: a = 1, b = -5, c = 4

  x0 = 0
  call solve_quad(a, b, c, x0, x, iters, .false.)
  print *, x, iters
  ! output: 1.0000000000000000, 5

  x0 = 7
  call solve_quad(a, b, c, x0, x, iters, .false.)
  print *, x, iters
  ! output: 4.0000000000000000, 6    
end
Rodrigo Rodrigues
  • 7,545
  • 1
  • 24
  • 36
  • 1
    I would have closed the question as a duplicate if it was just about this. We had pattern explained too many times already. Do you think this makes a difference for his second (and after the discussion the only real) question? If yes, it should be pointed out, because for the first simple problem using module variables instead of host association is perfectly valid as well. – Vladimir F Героям слава Sep 07 '18 at 07:01
  • I think this answer offers me another decent solution. I have never used the "internal procedure" so far. With it, I think I can simplify some of the source code without explicitly passing them. And I think this solution is different with the module method because each instance of the derived data type object is clearly separated. Thank you. Lesson learned. – Chang Sep 07 '18 at 21:48
  • The thing here is: this answer isn't different from the many others pointed by @VladimirF as possible duplicates; a careful reading of them would've lead you to a similar code. I just wrote this up because, despite that, this question didn't seem to get closed as duplicate and I don't want you or future readers of this Q to leave thinking there is no way of solving this. Summing up, pay close attention to advices from experienced ppl like him. – Rodrigo Rodrigues Sep 08 '18 at 02:21
1

After discussing with a colleague, I have a solution to my question 2.

If we have only one argument object for the solver module, then all the calling will access the same arguments because they share the same memory space.

To avoid this, we want to pass the argument object as an argument into the solver. So instead of using the default solver subroutine, we will re-write the Newton's method so it can accept additional argument.

(I used the simplest Newton subroutine earlier because I wanted to keep it untouched.)

In this way, we will define an array of argument objects and pass them during runtime.

Thank you for the comments.

Chang
  • 396
  • 4
  • 17