FORTRAN, подматрица как аргумент функции

kataich

Основной вопрос: что происходит, когда процедуре на вход подается подматрица?
Пример:

# cat main.f
IMPLICIT NONE
INTERFACE
SUBROUTINE PRIMER5(C)
REAL, DIMENSION(2:3,2:3)::C
END SUBROUTINE
END INTERFACE
REAL, DIMENSION(1:4,1:4)::D
D = 3
D(2:3,2:3) = -4
CALL PRIMER5(D(2:3,2:3
PRINT*, "D = ", D
END
# cat proc.f
SUBROUTINE PRIMER5(C)
IMPLICIT NONE
REAL, DIMENSION(2:3,2:3):: C

PRINT*,'C = ', C
C = 10
END
# gfortran main.f proc.f
# ./a.out
C = -4.0000000 -4.0000000 -4.0000000 -4.0000000
D = 3.0000000 3.0000000 3.0000000 3.0000000 3.0000000 10.000000 10.000000 3.0000000 3.0000000 10.000000 10.000000 3.0000000 3.0000000 3.0000000 3.0000000 3.0000000

Причина непоняток: в основной функции мы объявили массив D размерности 4x4 => в памяти он расположился непрерывным(?) куском по столбцам. Процедура PRIMER5 на вход тоже получает массив, который должен(?) располагаться непрерывным куском по столбцам. Когда мы вызываем CALL PRIMER5(D(2:3,2:3 мы передаем центральный квадратик 2x2, элементы которого в памяти не располагаются непрерывным образом. Однако, изменив значения массива C в процедуре, мы получили изменение и исходного массива. Если аргументы функции передаются по ссылке и передается указатель на элемент массива D, то тогда нарушается непрерывность массива в процедуре. Что происходит на самом деле?

Andbar

Операцию выделения подматрицы можно обратить.
SUBROUTINE PRIMER5(C)
IMPLICIT NONE
REAL, DIMENSION(2:3,2:3):: C
1:2,1:2 не пробовал?

spensnp

Если аргументы функции передаются по ссылке и передается указатель на элемент массива D, то тогда нарушается непрерывность массива в процедуре. Что происходит на самом деле?
почему вас это беспокоит?
это соответствует семантике языка, а остальное - дело компилятора, а не программиста фортраном.

kataich

1:2,1:2 не пробовал?

<strip>
SUBROUTINE PRIMER5(C)
IMPLICIT NONE
REAL, DIMENSION(1:2,1:2):: C
PRINT*, 'C(2,2) = ', C(2,2)
<strip>
CALL PRIMER5(D(2:3,2:3
<strip>

Все работает, как и ожидалось, точно также. Менять способ нумерации элементов массива вполне себе позволяется.

kataich

почему вас это беспокоит?

Потому что если происходит дополнительное выделение на стеке/в куче и т д, копирование туда-сюда, то я подумаю как сделать так, чтобы передавался один лишь указатель. Если же накладные расходы невелики, то буду использовать данный способ.
это соответствует семантике языка, а остальное - дело компилятора, а не программиста фортраном.
Если я не знаю (хотя бы примерно) механизмов, которые стоят за тем или иным синтаксическим выражением языка,
то я не уверен, что данное выражение делает то, что требуется мне. Тяжело каждый раз искать подтверждение в книге.
P.S. Если передается указатель

# cat main.f
IMPLICIT NONE
INTERFACE
SUBROUTINE PRIMER5(C)
REAL, DIMENSION(2:3,2:3)::C
END SUBROUTINE
END INTERFACE
REAL, DIMENSION(1:4,1:4)::D
D = 3
D(2:3,2:3) = -4
CALL PRIMER5(D(2,2))
PRINT*, "D = ", D
END

# cat proc.f
SUBROUTINE PRIMER5(C)
IMPLICIT NONE
REAL, DIMENSION(2:3,2:3):: C

PRINT*,'C = ', C
C = 10
END
# gfortran main.f proc.f
# ./a.out
C = -4.0000000 -4.0000000 3.0000000 3.0000000
D = 3.0000000 3.0000000 3.0000000 3.0000000 3.0000000 10.000000 10.000000 10.000000 10.000000 -4.0000000 -4.0000000 3.0000000 3.0000000 3.0000000 3.0000000 3.0000000

spensnp

P.S. Если передается указатель
ну это тоже ожидаемо )
ЗЫ
думаю что механизмы реализации сечений массивов запросто могут различаются в ifort, gfortran, g95 и наверняка там не обошлось без магических векторных инструкций и тайных знаниях о кэшировании данных процессорами... так что имхо лучше потестить производительность операций с сечениями на живых примерах и сделать вывод

yroslavasako

ни разу в lapack не видел передачу подматрицы (но я не все методы читал). Обычно передаётся кусок памяти (непрерывный и заранее выделенный если нужно индексы. Подматрицу - не fortran way

spensnp

> ни разу в lapack не видел передачу подматрицы
референсный lapack реализован с использованием стандарта fortran 77, в котором нет сечений
> Подматрицу - не fortran way[/quote]
сечения появились начиная с fortran 90

kataich

Спасибо за обсуждение.
Хотел бы задать еще один вопрос, но новый тред заводить не хочется.
Допустим есть две переменные:
INTEGER BOUND_LEFT, BOUND_RIGHT, обозначающие левую и правую границы массива. Значения их на момент компиляции не известны и определятся только на этапе исполнения. Однако хочется написать функцию,
на вход которой подается массив ARRAY( таким образом, чтобы нумерация массива внутри функции начиналась
бы с BOUND_LEFT и заканчивалась BOUND_RIGHT, то есть типа этого

SUBROUTINE FUNCTION(ARRAY)
REAL, DIMENSION(BOUND_LEFT: BOUND_RIGHT):: ARRAY
.............

Каков fortran way в этом случае? Заранее благодарю :)

kataich

Прошу покритиковать вот такой способ:

$ cat a.f
IMPLICIT NONE
INTERFACE
SUBROUTINE PRIMER5(C, M, N)
INTEGER M, N
REAL, DIMENSION(M:N)::C
END SUBROUTINE
END INTERFACE
REAL, DIMENSION(1:4)::D
INTEGER I, J
D = 3
D(2:3) = 1
I = 2
J = 3
CALL PRIMER5(D(2:3 I, J)
PRINT*, "D = ", D
END
$ cat b.f
SUBROUTINE PRIMER5(C, M, N)
IMPLICIT NONE
INTEGER M, N
REAL, DIMENSION(M:N):: C

PRINT*,'SIZE1 = ', SIZE(C, 1)
PRINT*,'C = ', C
C = 10
END


P.S. Я не волшебник - я только учусь, поэтому рад другим способам.

spensnp

Каков fortran way в этом случае?
я вовсе не гуру фортрана, и если уж на чистоту - я его вообще не перевариваю
зы пример выглядит логично, компилируется и работает, значит все ок )

slonishka

зы пример выглядит логично, компилируется и работает, значит все ок )
Why Language Features Matter. The Way the Fortran Bug Really Happened!
The details of a programming language really matter. They matter because the
details make the difference between a reliable language and an error-prone one.
This was dramatically revealed in Summer 1961 by a programmer at NASA, testing a
Fortran subroutine used to calculate orbital trajectories[1]. The subroutine had
already been used for several brief Mercury flights, but it was mysteriously not
providing the precision that was expected and needed for the forthcoming orbital
and lunar missions. The results were close, but not quite as accurate as
expected.
-1- The story is very widely misreported, and inaccurate versions appear in many
programming language texts. Indeed, it has become a classic urban legend
among programmers. The definitive account, from Fred Webb who worked at NASA
at the time and saw the actual source code, can be seen in "Fortran Story —
The Real Scoop" in Forum on Risks to the Public in Computers and Related
Systems, vol. 9, no. 54, ACM Committee on Computers and Public Policy,
December 12, 1989.
After checking the algorithm, the data, and the expected results at great
length, the engineer finally noticed this statement in the code:
DO 10 I=1.10
Clearly, the programmer had intended to write a `DO` loop of the form:
DO 10 I=1,10
In Fortran, blank characters are not significant and can even occur in the
middle of an identifier. The designers of Fortran intended this to help
cardpunch walloppers and to aid the readability of programs, so you could have
identifiers like `MAX Y`. Unfortunately, the compiler quite correctly read
the statement as `DO10I = 1.10`
Variables do not have to be declared in Fortran. The statement as written caused
the value 1.1 to be assigned to the implicitly declared floating point variable
`DO10I`. The statements in the body of the intended loop were executed once
instead of ten times, and a first approximation of a calculation took place
instead of an iterative convergence. After correcting the period to a comma, the
results were correct to the expected accuracy.
Оставить комментарий
Имя или ник:
Комментарий: