lat_radians = num[1:4000] (can also be df$lat_radians if easier) and solar_declination = num [1:365].

As a result, I want a data.frame that gives 'dl' for each lat_radians (row) multiplied by every solar_declination seperately (365 columns). So, iterations for full solar_declination[1:365]

The function itself works. However it only gives values for solar_declination[1].

Anyone who has a suggestion on how to accomplish this?

Every Rproblem can be though of with advantage as the interaction of three objects— an existing object, x , a desired object,y, and a function, f, that will return a value of y given x as an argument. In other words, school algebra— f(x) = y. Any of the objects can be composites.

Applied to this problem x is a composite object of two vectors lat_radians and solar_declination, the proto-y_0 is a lat_radians by solar_declinationmatrix object of dim [4000,365], one part of f is dl and the other part is a function to apply populate each element of y_0 by applying dl to each element to produce y.

Here's a toy example

my_func <- function(x,y) x*y
mat <- matrix(nrow = 3, ncol = 3, byrow = TRUE)
for (i in 1:3) {
for (j in 1:3) mat[i,j] = my_func(i,j)
}
mat
#> [,1] [,2] [,3]
#> [1,] 1 2 3
#> [2,] 2 4 6
#> [3,] 3 6 9