Yes, if you use NLsolve as in my post above with autodiff=:forward then it uses automatic differentiation to form the Jacobian matrix for Newton iterations. If you have a linear system of equations, the Jacobian matrix is exactly the matrix form of the problem, and Newton iterations will converge in a single \ step.
And if you know that your equations are linear, you can do the differentiation directly using the ForwardDiff.jl package:
julia> using ForwardDiff
julia> F(x) = [1- x[1] - x[2], 8 - x[1] - 3*x[2]]
F (generic function with 1 method)
julia> J = ForwardDiff.jacobian(F, [0,0])
2×2 Array{Int64,2}:
-1 -1
-1 -3
julia> x = -J \ F([0,0])
2-element Array{Float64,1}:
-2.5
3.5
Here, the Jacobian J corresponds to the “matrix form” of your problem for a right-hand-side of -F([0,0]), and appropriate use of \ (equivalent to a single Newton step) gives the same solution x as above.