2

How to verify in Fortran whether an iterative formula of a non-linear system will converge to the root near (x,y)?

It was easy for a programming language which support symbolic computations. But how to do that in Fortran? Like getting partial derivative of the component functions and check whether they bounded near the root. But I couldn't do that in fortran or haven't the idea how to do that. It will be a great help for me if anyone give me some idea for the following non-linear system now or if possible for a general case.

I want to use Fixed point iteration method for this case

Main system:

x^2+y=7
x-y^2=4

Iterative form (given):

X(n+1)=\sqrt(7-Y(n)),
Y(n+1)=-\sqrt(X(n)-4),
(x0,y0)=(4.4,1.0)

Theorem (which I follow)

image

The issue is, I need to check the boundedness of the partial derivatives of \sqrt(7-Y) and -\sqrt(X-4) on some region around (x0,y0)=(4.4,1.0). I can write the partial derivative function in fortran but how to evaluate so many values and check it is bounded around the (4.4,1.0).

Update

One possibly right solution would be to get arrays of values around (4.4,1.0) like (4.4-h,1.0-h)*(4.4+h,1.0+h) and evaluate the defined partial derivative function and approximate their boundedness. I haven't encounter such problem in Fortran, so any suggestion on that also can help me a lot.

Improve this question
4
  • What is your exact question? 1. How to "verify in Fortran" that some system converges. 2. How to verify, that some system will "converge in Fortran". 1 or 2? These are very different. Commented Oct 15, 2021 at 11:49
  • 1, @VladimirF. Even more specific, How to "verify in Fortran" that the given system converges. Commented Oct 15, 2021 at 11:51
  • Yes. I want to use Fixed point iteration method for this case @veryreverie Commented Oct 15, 2021 at 19:13
  • Is mentioning the theorem statement help? @veryreverie Commented Oct 15, 2021 at 20:13

2 Answers 2

Reset to default
0

If you just want to check the boundedness of a function on a grid, you can do something like

program verify_convergence
  implicit none
  
  integer, parameter :: dp = selected_real_kind(15, 307)
  
  real(dp) :: centre_point(2)
  real(dp) :: step_size(2)
  integer :: no_steps(2)
  real(dp) :: point(2)
  real(dp) :: derivative(2)
  real(dp) :: threshold
  
  integer :: i,j
  real(dp) :: x,y
  
  ! Set fixed parameters
  centre_point = [4.4_dp, 1.0_dp]
  step_size = [0.1_dp, 0.1_dp]
  no_steps = [10, 10]
  threshold = 0.1_dp
  
  ! Loop over a 2-D grid of points around the centre point
  do i=-no_steps(1),no_steps(1)
    do j=-no_steps(2),no_steps(2)
      ! Generate the point, calculate the derivative at that point,
      !    and stop with a message if the derivative is not bounded.
      point = centre_point + [i*step_size(1), j*step_size(2)]
      derivative = calculate_derivative(point)
      if (any(abs(derivative)>threshold)) then
        write(*,*) 'Derivative not bounded'
        stop
      endif
    enddo
  enddo
  
  write(*,*) 'Derivative bounded'
contains
  
  ! Takes a co-ordinate, and returns the derivative.
  ! Replace this with whatever function you actually need.
  function calculate_derivative(point) result(output)
    real(dp), intent(in) :: point(2)
    real(dp) :: output(2)
    
    output = [sqrt(7-point(2)), -sqrt(point(1)-4)]
  end function
end program

I know the function calculate_derivative doesn't do what you want it to, but I'm not sure what function you actually want from your question. Just replace this function as required.

Improve this answer
Sign up to request clarification or add additional context in comments.

1 Comment

It gives me some idea @veryreverie. Let me test this on my machine. Thanks for your answer.
-1

The main question is different: How can you calculate the solution of the mathematical problem without the help of any software? If you know that, we can program it in fortran or any language.

In particular, and assuming that n=0,1,2,3... to solve your problem you need to know X(0) and Y(0). With this you calculate

X(1)=\sqrt(7-Y(0)), Y(1)=-\sqrt(X(0)-4)

Now you know X(1) and Y(1), then you can calculate

X(2)=\sqrt(7-Y(1)), Y(2)=-\sqrt(X(1)-4)

etc.

If your system of equations converge to something, until some iterations (for example, n=10, 20, 100) you going the check that. But because the nature of fortran, it will not give you the solution in a symbolic way, it is not its goal.

Improve this answer

2 Comments

(x0,y0) is given as (4.4,1.0). I am not asking how to code to get iteration value. I was asking how to check the system would converge near that point @ramrebol
I insist that to know that the algorithm is correct (it converges) you must verify it with pencil, paper and patience. What you can do in fortran is code it and see how it behaves for n = 10, n = 20, n = 100 ... and if it seems to converge as n is larger, then you are verifying that it is correct. But, of course, it is not a demonstration. Fortran (formula translator) allows you to give orders for the computer to calculate, it does not make any demonstration.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.