The shape of the heteroscedasticity in the scatterplot is a strong indicator that a nonlinear re-expression of $y$, the school grade point total, will cure the problem.
A good choice at the outset would be the "folded logarithm" or, equivalently, a logit-like formula
$$z = \log(y / (100 + \varepsilon - y))$$
where $\varepsilon$ is a small positive value, perhaps around $1/2$ or $1$, chosen to avoid singular calculations.
To demonstrate this, I simulated data like the ones shown in the question. The x-coordinates appear concentrated at the right side, so (in R) I created 500 values by scaling a Beta$(4,1)$ distribution to the range $[0,100]$:
set.seed(17)
x <- 100 * rbeta(500, 4, 1)
Next I posited that the logits of the school scores, $z$, would have a linear relationship to the poverty indexes $x$ plus normally-distributed random error. I chose an intercept and slope to reproduce the apparent fit and a standard deviation of the error to reproduce the scatter:
z <- 4 - .025*x + rnorm(length(x), sd=1)
logistic <- function(y) 1 / (1 + exp(-y))
y <- 100 * logistic(z)
Let's look at the plot to see how successful these formulae are at emulating the data:
par(mfrow=c(1,1))
plot(x,y, main="Simulated Data")
That looks quite similar to the plot included in the question.
Now let's display the re-expressed scores:
logit <- function(x) log(x/(100-x))
plot(x, logit(y), ylab="Logit(y)", main="Re-expressed Data")
lines(lowess(logit(y) ~ x), col="Red", lwd=2)
To the uninitiated, this scatterplot may appear heteroscedastic, but that is an illusion brought about by the concentration of $x$ values at the right: as is clear from the code, these data are perfectly homoscedastic. To make this clear, construct a wandering schematic plot by slicing $x$ at key quantiles and displaying side-by-side boxplots of the residuals:
thresh <- c(0, 1/16, 1/8, 1/4, 1/2, 3/4, 7/8, 15/16, 1)
boxplot(residuals(lm(logit(y) ~ x)) ~ cut(x, quantile(x, probs=thresh)),
main="Residuals of logit(y)~x")
These plots show the residuals have essentially the same dispersion across the entire range of $x$ values, demonstrating the homoscedasticity of the re-expressed scores.
Finally--to drive these points home--let's compare least-squares fits of the original scores and of their re-expressed values. First, the original scores:
summary(lm(y ~ x))
Estimate Std. Error t value Pr(>|t|)
(Intercept) 101.8218 2.6770 38.04 < 2e-16 ***
x -0.2154 0.0327 -6.58 1.2e-10 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 12 on 498 degrees of freedom
Multiple R-squared: 0.0801, Adjusted R-squared: 0.0782
F-statistic: 43.4 on 1 and 498 DF, p-value: 1.16e-10
Standard diagnostic plots (not shown here) are terrible, showing strongly heteroscedastic residuals, non-normal residuals, etc.
Now for the re-expressed data:
summary(lm(logit(y) ~ x))
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3.75714 0.21751 17.27 < 2e-16 ***
x -0.02157 0.00266 -8.12 3.8e-15 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.975 on 498 degrees of freedom
Multiple R-squared: 0.117, Adjusted R-squared: 0.115
F-statistic: 65.9 on 1 and 498 DF, p-value: 3.76e-15
The original coefficients of $4$ and $-0.025$ are adequately estimated as $3.75$ and $-0.022$, respectively. The residual standard error of $0.975$ is close to the intended SD of $1$. All diagnostics are beautiful.
I am looking at state-wide data (entire population) of a school's grade as a function of the school's poverty index. The data appears to me to be an unconditional heteroskedastic distribution. I am drawing a line of best fit (LBF) through the data using linear regression (I changed this from 2-degree poly). What I am attempting to do is look at how each school is doing compared to how it was ‘predicted’ to do. In past years, a different test was used and the result was that the data was clearly linear, and moderately heteroskedastic. I then used SD as a measure so that I could look across multiple years. 
