I'm not sure I understand your problem. What exactly did you expect your result will be with this code? Note that x is not sorted, and you're plotting the points against 1:10000.
Here's an way to plot with ggplot2. I simulated in the range [-5, 5], but you can change as you prefer.
Note that the implementations in the packages fGarch and sn vary.
library(fGarch)
#> Loading required package: timeDate
#> Loading required package: timeSeries
#> Loading required package: fBasics
library(ggplot2)
library(sn)
#> Loading required package: stats4
#>
#> Attaching package: 'sn'
#> The following object is masked from 'package:fBasics':
#>
#> vech
#> The following object is masked from 'package:stats':
#>
#> sd
central_tendency_sim <- 0
dispersion_sim <- 1
skewness_sim <- 1.5
N_sim <- 10000
obs_sim <- seq(from = -5,
to = 5,
length.out = N_sim)
ggplot(data = data.frame(u = obs_sim),
mapping = aes(x = u)) +
stat_function(mapping = aes(colour = "fGarch"),
fun = dsnorm,
args = list(mean = central_tendency_sim,
sd = dispersion_sim,
xi = skewness_sim)) +
stat_function(mapping = aes(colour = "sn"),
fun = dsn,
args = list(xi = central_tendency_sim,
omega = dispersion_sim,
alpha = skewness_sim)) +
stat_function(mapping = aes(colour = "standard"),
fun = dnorm) +
scale_colour_manual(name = NULL,
values = c("red", "blue", "green"))

Created on 2019-09-06 by the reprex package (v0.3.0)
In case you want to plot in base R, that also works.
library(fGarch)
#> Loading required package: timeDate
#> Loading required package: timeSeries
#> Loading required package: fBasics
library(sn)
#> Loading required package: stats4
#>
#> Attaching package: 'sn'
#> The following object is masked from 'package:fBasics':
#>
#> vech
#> The following object is masked from 'package:stats':
#>
#> sd
central_tendency_sim <- 0
dispersion_sim <- 1
skewness_sim <- 1.5
N_sim <- 10000
obs_sim <- seq(from = -5,
to = 5,
length.out = N_sim)
den_sim_fGarch <- dsnorm(x = obs_sim,
mean = central_tendency_sim,
sd = dispersion_sim,
xi = skewness_sim)
den_sim_sn <- dsn(x = obs_sim,
xi = central_tendency_sim,
omega = dispersion_sim,
alpha = skewness_sim)
par(mfrow = c(2, 2))
plot(x = obs_sim,
y = den_sim_fGarch,
type = "l",
xlab = "observations",
ylab = "densities",
main = "plot points and join - fGarch",
col = "violetred")
curve(expr = dnorm,
from = -5,
to = 5,
n = N_sim,
add = TRUE,
col = "limegreen",
lty = 2)
plot(x = obs_sim,
y = den_sim_sn,
type = "l",
xlab = "observations",
ylab = "densities",
main = "plot points and join - sn",
col = "navyblue")
curve(expr = dnorm,
from = -5,
to = 5,
n = N_sim,
add = TRUE,
col = "limegreen",
lty = 2)
plot.function(x = function(t) dsnorm(x = t,
mean = central_tendency_sim,
sd = dispersion_sim,
xi = skewness_sim),
to = 5,
from = -5,
n = N_sim,
xlab = "observations",
ylab = "densities",
main = "plot function - fGarch",
col = "palevioletred")
curve(expr = dnorm,
from = -5,
to = 5,
n = N_sim,
add = TRUE,
col = "limegreen",
lty = 2)
plot.function(x = function(t) dsn(x = t,
xi = central_tendency_sim,
omega = dispersion_sim,
alpha = skewness_sim),
to = 5,
from = -5,
n = N_sim,
xlab = "observations",
ylab = "densities",
main = "plot function - sn",
col = "skyblue")
curve(expr = dnorm,
from = -5,
to = 5,
n = N_sim,
add = TRUE,
col = "limegreen",
lty = 2)

dev.off()
#> null device
#> 1
Created on 2019-09-06 by the reprex package (v0.3.0)